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ABSTRACT 


The  entraining-parcel  model  of  the  rise  and  expansion  of  the  nuclear 
cloud  is  revised  to  allow  for  (l)  production  of  turbulent  kinetic  energy 
from  kinetic  energy  of  rise  during  the  momentum-conserving,  inelastic- 
collision  entrainment,  so  that  total  kinetic  energy  is  conserved  in  en¬ 
trainment.  This  production  of  turbulent  energy  is  in  addition  to  that 
due  to  eddy  viscosity;  (2)  dissipation  of  turbulent  energy  to  heat. 

The  resulting  nuclear- cloud  model  is  represented  as  an  energy 
cycle  between  enthalpy  and  energy  of  rise,  turbulent  energy,  and  poten¬ 
tial  energy. 

Calculations  of  toroidal  circulation  and  particle  notion  in  the 
nuclear  cloud  are  discussed,  using  the  cloud  model  and  turbulent  similarity 
theory.  If  the  toroidal  circulation  is  to  be  represented  by  a  vertex 
ring  superimposed  on  the  parcel-model  cloud,  the  vortex  ring  can  be  con¬ 
sidered  the  largest  eddy  in  the  turbulent  spectrum,  containing  a  fixed 
fraction  of  the  turbulent  energy.  This  gives  an  estimate  of  the  circu¬ 
lation  without  recourse  to,  but  in  general  agreement  with,  estimates 
based  on  measurements  of  nuclear  cloud  films.  Alternately,  circulation 
may  be  calculated  by  an  adaptation  of  Kelvin's  theorem.  This  method 
assumes  neither  a  particular  vortex  form  nor  steady- state  flow,  and  is 
therefore  more  consistent  with  actual  cloud  conditions. 

Attempts  have  been  made  to  calculate  the  effect  of  toroidal  circu¬ 
lation  on  particle  dispersion  from  the  cloud  using  laminar  flow  methods. 
Here,  this  dispersion  is,  instead,  represented  as  due  to  turbulent  dif¬ 
fusion,  using  the  calculated  dissipation  rate  as  the  governing  parameter. 
Diffusion  coefficients  and  concentration  gradients  are  derived  from  tur¬ 
bulent  similarity  theory.  The  resulting  dispersion  rate  is  shown  to  be 
small  compared  with  that  due  to  gravitational  fallout  rate  indicating 
that  dispersion  induced  by  circulation  or  turbulence  can  be  ignored. 


SUMMARY  PAGE 


The  Problem 


A  vortex  ring  or  so-called  toroidal  circulation  is  often  observed 
in  the  intensely  turbulent  cloud  rising  from  a  nuclear  explosion.  Par¬ 
ticles  of  different  sizes  in  this  nuclear  cloud  may  follow  the  circula-° 
tory  motion  of  the  fluid  to  different  extents,  and  some  of  them  may  be 
centrifuged  out  of  the  cloud.  The  circulation,  then,  affects  the  dis¬ 
persion  of  particles  from  the  cloud  and  the  distribution  of  radioactive 
fallout  in  the  air  and  on  the  ground. 

Calculation  of  the  fluid  motion  and  the  particle  trajectories  by 
classical  hydrodynamic  methods  presents  both  theoretical  and  computational 
difficulties,  due  largely  to  the  turbulent  character  of  the  flow  in  and 
around  the  cloud.  Calculations  of  cloud  motion  can  be  considerably 
simplified  by  using  a  parcel  method,  which  gives  average  values  of  cloud 
velocity,  size  t\nd  other  variables  as  functions  of  time.  An  entralning- 
parcel  cloud  model  has  previously  been  developed  at  this  laboratory. 

The  objectives  of  the  present  work  were: 

(1)  to  revise  the  cloud  model  to  take  into  account  dissipation 

of  turbulent  energy  to  heat,  and  conservation  of  kinetic  plus  turbulent 
energy  in  entrainment. 

(2)  to  discuss  attempts  to  calculate  toroidal  circulation  in  the 
nuclear  cloud  and  its  effects  on  particle  dispersion. 

(3)  to  investigate  alternate  methods  of  modeling  non-gravitational 
particle  dispersion,  using  turbulence  theory. 


Findings 

(1)  In  the  revised  cloud  model,  total  energy  is  conserved.  The 
model  implies  a  nuclear  cloud  energy  cycle. 

(2)  If  toroidal  circulation  is  to  be  represented  by  a  classical 
vortex  ring  superimposed  on  a  parcel-model  cloud,  then  the  vortex  para¬ 
meters  may  be  estimated  from  the  cloud  model  consistent  with  turbulenc 
theory.  The  calculated  values  of  circulation  are  in  general  agreement 
with  the  few,  scattered  values  based  on  measurements  of  nuclear  cloud 
films.  Alternately,  the  rate  of  change  of  circulation  may  be  calculated 
by  a  free  adaptation  of  Kelvin's  theorem,  and  the  rate  integrated  to 
give  circulation.  This  method  does  not  require  assuming  a  particular 
vortex  structure  and  steady- state  flow.  It  is  therefore  more  consistent 
with  physical  conditions  in  the  nuclear  cloud. 

(3)  Non-gravitational  dispersion  of  particles  can  be  represented 
as  due  to  turbulent  diffusion  instead  of  laminar  vortex  flow.  Using 
turbulent  diffusion  coefficients  derived  from  small- scale- turbulence 
theory,  and  the  turbulent  energy  dissipation  rate  given  by  the  cloud 
model,  it  is  shown  that  turbulent- diffusive  dispersion  of  particles 

is  unimportant  compared  to  gravitational  fallout. 
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1.  INTRODUCTION 


1.1  Background 

One  of  the  striking  features  of  nuclear  explosions  is  the  vortex 
ring  and/or  "toroidal  circulation"  in  the  rising,  intensely  turbulent 
nuclear  clouf,.  It  has  been  suggested  that  particles  of  different  sizes 
in  the  cloud  follow  the  motion  of  the  fluid  circulation  to  different 
extents,  so  that  the  circulation  influences  the  dispersion  of  particles 
from  the  cloud  (b>  a  sort  of  centrifuging  process)  and  thus  the  distribu¬ 
tion  of  radioactive  fallout  in  the  air  and  on  the  ground.  Calculation 
of  the  fluid  motion  and  of  the  particle  trajectories  by  classical  hydro¬ 
dynamic  methods  presents  both  theoretical  and  c crnputat ional  difficulties, 
not  yet  overcome.  These  difficulties  arise  largely  from  the  turbulent 
character  of  the  flow  in  and  around  the  cloud. 

1.2  Objectives 

The  purposes  of  this  report  are 

1.  To  present  a  revision  of  a  previously  published1  entraining- 
parcel  model  of  the  atomic  cloud,  allowing,  as  before,  for  transfer  of 
kinetic  energy  of  cloud  rise  to  turbulence  by  eddy-viscosity,  and  alBO 
for 

(a)  additional  transfer  of  kinetic  energy  of  rise  to  turbulence 
by  the  inelastic-collision,  momentum- conserving  entrainment  process,  so 
that  total  energy  is  conserved 

(b)  ti-ansfer  of  energy  from  turbulence  to  heat. 

2.  To  discuss  the  energy  cycle  of  the  atomic  cloud  which  is  im¬ 
plicit  in  this  model. 

3.  To  discuss  attempts  to  superimpose  a  classical  vortex  motion  on 
a  par  cel-  method  ■  cloud,  and  to  show  how  such  attempts  can  at  least  be  made 
consistent  with  the  energy  balance  of  tbe  cloud,  and  hcrw  cloud  circula¬ 
tion  can  be  estimated  without  use  of  the  superimposed  vortex. 
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4.  To  show  how  the  particle-dispersive  effect  attributed  to  fluid 
circulation  can  be  represented  as  due  to  turbulent  diffusion  consistent 
with  the  parcel  method,  and  to  estimate  the  relative  importance  of  dif¬ 
fusive  and  gravitational  dispersion  of  particles  from  the  cloud. 

1 . 3  Approach 

The  principal  tool  used  in  the  present  study  is  the  revised  entraining- 
parcel  model  of  the  nuclear  cloud.  This  model  was  originally  developed 
for  water- surface  bursts,  but  is  also  applicable  to  air  bursts,  since, 
in  both  cases,  the  mass  fraction  of  condensed  matter  is  small.  The 
model  can  further  be  applied  to  land- surface  bursts  provided  the  mass 
and  size  distribution  of  condensed  matter  at  start,  of  cloud  rise  are 
specified. 

The  present  development  of  the  original  mo<?’l  .U  largely  inspired 
by  tir  bulence  theory,  especially  by  the  concept  of  uhe  energy  cascade 
and  the  theory  of  local  similarity,  or  similarity  of  small  eddies. 

Therefore,  the  relevant  parts  of  turbulence  theory  are  reviewed  in  this 
report.  (Sec.  3) 


2.  MODELS  OF  THE  NUCLEAR  CLOUD 

2.1  The  Formation  of  the  Nuclear  Cloud 

In  a  nuclear  erploaion,  typically,  from  1  KT  to  10  MT  of  energy 
(1012  to  lO1^  caloriea)  is  released  in  less  then  1  microsecond  (p  sec.) 
in  a  space  1  meter  in  diameter. 

Since  the  energy  release  ye r  unit  mass  is  thousands  of  times  that 
of  a  chemical  explosion,  the  resulting  "fireball"  consisting  of  vaporized 

7 

bomb  materials  is  initially  at  a  temperature  on  the  order  of  10  K  instead 
of  3000  -  5000  K  as  for  a  chemical  explosion,  and  a  pressure  of  several 
million  atmospheres. 

The  total  energy,  W,  of  a  nuclear  explosion  in  the  air,  near  sea 
level,  is  partitioned  thus: 

5$  =  initial  nuclear  radiation  (radiation  during  the  first 
minute  after  the  explosion) 

10$  =  residual  nuclear  radiation  (radiation  after  the  first 
minute  after  the  explosion) 

50$  =  blast  and/or  shock 
35$  =  thermal  radiation 

All  of  this  energy  is  eventually  degraded  to  heat. 

The  diameter  of  the  fireball  increases  by  radiative  and  shock 
heating  of  the  surrounding  air,  vhile  the  pressure  in  the  fireball  drops 
by  expansion  and  radiative  cooling.  Within  a  few  seconds  (about  1 
second  for  W  «  20  KT)  pressure  has  dropped  to  one  atmosphere  and  temper¬ 
ature  han  dropped  to  about  3000  to  4500  K.  The  fireball  then  starts 
rising  and  can  be  called  the  nuclear  cloud. 

At  start  of  cloud  rise,  the  fireball  contains  about  ^  of  the  total 
energy  of  the  explosion.  This  ^  value  need  not  be  the  same  as  the 
35$  thermal  energy  mentioned  above,  since  seme  of  Tthe  is  shock 
heating  and  seme  of  the  35$  vas  radiated  beyond  the  fireball. 

The  fraction,  f,  of  the  total,  energy  contained  in  the  fireball  at 
start  of  rise,  fW,  increases  slightly  irith  yield  but  the  subsequent  cloud 

hiBtory  is  not  sensitive  to  variations  in  assumed  initial  cloud  energy, 

4  i 

fit  .  Taylor  calculated  f  >  +  as  the  basis  of  shock  heating  plus  adiabatic 


expansion  to  atmospheric  pressure,  and  suggested  that  the  effects  of 
neglecting  both  energy  loss  by  radiation  and  the  increased  specific 
heat  of  air  at  high  temperature  tend  to  cancel  each  other. 

The  maximum  radius  of  the  fireball  before  it  starts  to  rise  and 
becomes  the  atomic  cloud  is  proport ional^  to  w2/5.  The  exponent  2/5 
is  close  to  l/3  and  the  fireball  size  is  approximately  as  if  1/3  of  the 
total  explosion  energy  were  used  to  heat  air  at  constant  pressure  to 
3000K.  The  calculated  fireball  volume  (and  thus  radius)  at  start  of 
rise  is  nearly  independent  of  assumed  fireball  temperature. 

(Proof:  using  the  perfect  gas  law,  and  the  equality  of  cloud  and 
ambient  pressure 

p  ■  pRT  =  p  RT  ;  fW  *>  me  (T-T  )  so  that 
ee  p  e 

„  m  R  fWT  _  fRW 

V  "J"?  c'fl-T  )  -7J 
p  e  p 

Here  the  subscript  e  denotes  initial  environment  conditions,  p  is  pres¬ 
sure,  p  is  density,  R  is  the  universal  gas  constant,  T  is  temperature, 

V  is  cloud  volume,  m  is  cloud  mass,  and  c^  is  specific  heat  at  constant 
pressure. ) 

At  start  of  rise,  the  cloud  may  contain  some  refractory  particulate 
matter;  more  and  larger  particles  form  in  the  rising,  cooling  cloud  by 
condensation  and  coagulation.  The  dispersion  of  these  particles  from 
the  cloud  may  be  affected  by  conditions  inside  the  cloud,  including  vor¬ 
tex  flow  or  turbulence.  The  following  section  considers  calculation  of 
such  conditions;  l.e.  cloud  models. 

2.2  Cloud  Calculations 

Calculations  of  the  behavior  of  tho  nuclear  cloud  have  been  made  by 
a  number  of  investigators  using  different  approaches,  often  patterned 
after  studies  of  cumulus  clouds  or  of  the  rise  of  hot  gases  from  fires 
and  factory  chimneys. 
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Possible  models  of  cloud  behavior  can  be  divided  into  those  using 
(l)  local  methods  and  (2)  parcel  methods. 


2.2.1  Local  Methods 

In  a  "local"  method,  the  flow  conditions  in  the  cloud  and  surrounding 
atmosphere  are  described  by  a  set  of  partial  differential  equations 
(Navier- Stokes  equations)  giving  local  values  of  velocity,  temperature, 
etc.,  as  functions  of  time.  For  numerical  computation,  the  partial 
differential  equations  are  replaced  by  a  set  of  ordinary  differential 
equations  at  each  point  of  a  grid  or  mesh  of  points  covering  the  flow 

e  f, 

region.  Some  beginnings  using  this  treatment  have  been  made^  for 
developing  cumulus  clouds,  making  the  assumption  that  temperature  diff¬ 
erences  between  cloud  and  atmosphere  are  small  relative  to  absolute 
temperature,  and  using  an  eddy  viscosity  in  the  Navier-Stokes  equations. 
The  values  of  velocity,  etc.,  computed  represent  time  averages,  not  the 
instantaneous  values  which  are  subject  to  turbulent  fluctuations. 

Eddy- viscosity  is  needed  in  the  equations  used  for  computation 
because  molecular  viscosity  (ordinary  viscosity)  operates  to  convert 
the  energy  of  only  the  smallest  eddies  of  a  turbulent  flow  into  heat. 

It  can  be  shown  (see  for  instance,  Section  3.1.5)  that  the  ratio  of  the 
size  of  the  smallest  eddies,  Xq,  to  the  overall  flow  dimensions,  &,  for 
one-di-mens ional  flow,  is 


where  Re  is  the  Reynolds  number  of  the  flow.  Such  a  flow  could  be 

pictured  as  containing  a  rcw  of  l/\  such  eddies,  (on  which  are  super- 

°  3/4 

imposed  larger  eddies)  and  so  having  Re  '  degrees  of  freedom.  Then 

3/2 

the  mesh  for  a  two-dimensional  calculation  must  have  Re  '  points,  so 
that  numerical  computation  is  out  of  the  question  when  Re  is  large.  The 
eddy  viscosity  is  used  to  blur  together  small  eddies,  and  thus  reduce 
the  number  of  degrees  of  freedom  of  the  flow  problem,  and  the  required 
number  of  mesh  roints,  to  where  numerical  computation  is  possible. 


5 


No  local- method  computation  has  yet  been  carried  out  for  nuclear 
clouds.  Even  for  cumulus  clouds,  the  original  two-dimensional  array  of 
mesh  points  (in  an  axisymmetric  finite- difference  calculation)  is  rapidly 
distorted.  For  nuclear  clouds,  this  distortion  would  be  accompanied  by 
extreme  temperature  gradients,  so  that  heat  conduction  could  not  be 
neglected.  Presumably,  ar.  eddy  conductivity  would  be  required  as  well 
as  eddy  viscosity.  Altogether,  one  is  inclined  to  question  the  use  of 
eddy  transport  coefficients  to  calculate  mean  local  conuitions  at  the 
extremely  high  Reynolds  numbers  (~10^)  of  these  intensely  turbulent 
flows.  In  any  case,  the  attempt  has  not  yet  been  made. 

2.2.2  Farcel  Methods 

In  a  parcel  method,  the  cloud  is  treated  as  a  whole,  as  if  all 
parts  of  it  had  the  same  properties.  Temperature,  velocity,  etc.,  are 
represented  by  average  values  for  the  whole  cloud.  Pressure  is  taken  as 
equal  to  ambient  pressure  at  the  altitude  of  the  center  of  the  cloud. 

Thus  a  parcel  method  is  a  simple,  cheap  substitute  for  a  local  method. 

But  it  is  far  more  amenable  to  numerical  computation.  And  with  a  few 
enlightened  assumptions,  a  parcel  method  can  give  a  wide  variety  of 
information,  including  an  estimate  of  particle  dispersion  by  toroidal 
circulation,  represented  es  turbulent  diffusion.  See  Sec.  5. 

2.3  The  TR-?4l  Atomic  Cloud  Model,  Revised 

An  atomic  cloud  model,  using  an  entraining  parcel  method  has  been 
developed  at  this  Laboratory1.  This  model,  (referred  to  as  the  TR-7^1 
model)  was  particularly  designed  for  sea-water-surface  nuclear  explosions. 
The  model  treats  the  cloud  gas  as  a  mixture  of  air  and  water  va>.,or  and 
allows  for  the  release  of  latent  heat  by  water-vapor  condensation.  The 
key  assumptions  in  this  model  are: 
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1.  The  effective  rate  of  flow  of  ambient  air  into  the  cloud  (en¬ 
trainment)  per  unit  surface  area  is  equal  to  the  product  of  a 
constant,  \,  a  characteristic  velocity,  v,  and  the  ratio  of 
cloud  density  to  ambient  density,  p/pg*  Some  previous  cloud 
studies  have  used  the  same  entrainment  rule  except  for  omitting 
the  density  ratio. 

2.  Cloud  rise  is  retarded  by  both  entrainment  and  an  apparent  eddy- 
viscous  force,  directly  proportional  to  the  same  characteristic 
velocity  and  inversely  proportional  to  the  density  ratio.  (This 
is  formally  equivalent  to  a  drag  coefficient.) 

3.  The  kinetic  energy  of  rise  corresponding  to  the  momentum  lost 

by  eddy  viscosity  is  converted  into  kinetic  energy  of  turbulence, 
which  remains  in  the  cloud. 

4.  The  characteristic  velocity,  v,  is  the  greater  of  absolute  rate 
of  cloud  rise,  u,  and  average  velocity  of  turbulence,  VST*. 

Thus  entrainment  does  not  necessarily  end  when  cloud  rise  ends, 
but  continues  as  a  turbulent  diffusion. 

5.  In  accelerating  from  rest  to  its  maximum  velocity,  the  cloud 
must  set  in  motion  a  volume  of  ambient  air  equal  to  one  half 
the  initial  cloud  volume.  Because  of  this  "virtual  mass",  the 
acceleration  is  always  less  than  twice  that  of  gravity. 

No  provision  was  made  in  this  model  for  transformation  of  turbulent 

energy  into  heat.  Not  only  does  such  transformation  actually  take  place, 

but  also  it  may  affect  the  formation  of  fallout  particles  through  the 

7 

mechanism  of  turbulent  coagulation.  The  rate  of  this  transformation, 
the  so-called  "dissipation  rate",  £  ,  is  found  both  theoretically  and 
experimentally  to  be  proportional  to  the  cube  of  a  large-scale  turbulent 
velocity  divided  by  a  large-scale  length,  l.  (see  Section  3>  Turbulence) 
Therefore,  we  now  add  one  assumption  to  the  cloud  model. 

6.  Turbulent  energy  is  dissipated  at  a  rate  proportional  to  the 
cube  of  the  average  velocity  of  turbulence  divided  by  a  ch&r- 

*  The  symbol  was  used  in  TR-7^1  for  the  quantity  called  E  in  the 
present  report. 


7 


acteristic  length.  (See  Revisions,  Appendix  A)  The  average 
velocity  of  turbulence  is  defined  as  yjz E  where  E  is  the  tur¬ 
bulent  energy  per  unit  mass.  'The  characteristic  length  used 
is  the  vertical  radius  of  the  cloud. 

It  was  not  previously  noted  that  in  entrainment  of  stationary  am¬ 
bient  air  with  conservation  of  momentum,  kinetic  energy  is  lost,  i.e., 
entrainment  is  an  inelastic-collision  process.  The  proof  s  as  follows: 
The  rate  of  change  of  cloud  momentum  per  unit  mass,  u, 

(i.e.,  velocity)  with  time,  t,  due  to  entrainment,  is 

du  1  dm  _  . , 


du  1  dm 

dt  U  m  dt 


(2.3.1) 


The  corresponding  change  in  kinetic  energy  per  unit  mass  is 

du  d(u2/g)  =  .  pf— )  i  —  (p?2) 

dt  dt  \2  I  m  dt  W-3'*) 

On  the  other  hand,  if  total  kinetic  energy  is  conserved,  then 
for  a  small  change  in  velocity  and  mass 
12  1  2 

—  mu  =  —  (m+-dm)  (u+du) f  and  neglecting  2nd-order  terms 


d(u2/2)  _  u2  1  dm 

dt  2  m  dt 


(2.3.3) 


Entrainment,  then,  Results  in  a  loss  of  total  kinetic  erergy,  at  a 
rate,  per  unit  mass,  of  j-j  The  law  of  conservation  of  total  energy 

and  the  concept  of  turbulent  energy  already  used  suggest  the  following 
assumption : 

7 .  The  kinetic  energy  of  rise  lost  in  the  inelastic- 
collision,  momentum- conserving  entrainment  process, 
remains  in  the  cloud  as  'turbulent  energy. 

A  simplified  version  of  the  essential  equations  is  given  below, 
neglecting  the  effects  of  the  2-gas  mixture,  latent  heat  release,  initial 
virtual  mass  and  condensed-water  mass  fraction  of  the  cloud.  These  effects 
sure  considered  in  the  full  set  of  equations  used  in  computation.  (See  TR- 
741  and  Appendix  A) 

Massj.  if£  =  |^v  (2.3.4) 
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where 


m  =  mass 


t  =  time 

S  =  cloud  surface  area 
V  =  cloud  volume 

\  =  dimensionless  entrainment  constant 
v  =  characteristic  velocity  of  entrainment 


where  u  =  rate  of  rise 

T  =  cloud  temperature 
Tg=  environment  temperature 
g  =  acceleration  of  gravity 

dimensionless  drag  or  eddy  viscosity  constant 
/  =  a  characteristic  length  of  the  cloud 


Momentum:  ~  =  ^ - 1^ 


=  _  fL 

dt  c 


(T-T  )  ±  ^ 
e  m  dt 


1  dm  c 


(2.3.6) 


where 


specific  heat  of  air 


Turbulent  kinetic  energy  density: 


dE  T  v  2 

dt  =  ^2  T  J  u 

6 


+  u_ldm  ldm  (?  ■>  n\ 

+  T  m  dt  E  m  dt  E  (2'3*7) 


where  E  =  turbulent  energy  per  unit  mass 
u2  1  dm 

The  term  —  —  —  compensates  for  the  loss  of  kinetic  energy  of  rise 
d  m  dt 

in  entrainment  (assumption  7)*  Thus,  some  turbulent  energy  is  produced 
even  in  the  absence  of  eddy  viscosity,  i.e.  even  if  kg  =  0. 

Environment  temperature,  T  ,  is  a  specified  function  of  height,  z, 
which  is  given  by 

(2.3.8) 
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Dissipation  rate  is  given  by 

e  =  k,  (2E)^2/f  (2.3.9) 

.) 

where  is  a  dimensionless  constant. 

Cloud  form  is  a  snhere  (initially  tangent  to  sea  level,  say,  in  the 
case  of  a  surface  burst)  until  the  top  (not  the  center,  as  in  TR-74l) 
reaches  the  tronopause,  and  is  a  horizontally  expanding  spheroid  thereafter. 
The  vertical  radius  of  the  spheroid  is  fixed  as  the  sphere's  radius  at 
the  tropopause.  Hie  cloud  volume  and  surface  area  can  then  be  calculated 
from  the  perfect  gas  law  using  the  assumption  of  pressure  equilibrium. 

The  model  actually  uses  a  differential  equation  for  volume,  so  that  the 
gas  law  is  available  as  a  cross-check. 

The  characteristic  length,  i  ,  is  taken  as  the  (vertical)  radius  of 
the  cloud,  'lhe  characteristic  velocity,  v,  is  taken  as  max  (|u|,'J2E), 

The  equations  have  been  programmed  for  machine  computation,  and  the 
computations  give  predictions  for  rate  of  rise,  cloud  size,  final  cloud 
height  and  the  late  horizontal  expansion  of  high-yield  clouds,  in  general 
agreement  with  observations  of  nuclear  clouds.  Late  horizontal  expansion 
had  not  previously  been  predicted  by  a  cloud  model.  Detailed  revisions 
to  the  TR-741  cloud  model  are  given  in  Appendix  A.  Appendix  F  gives  the 
corresponding  revised  computer  program.  Tables  2.1  and  2.2  are  sample 
computer  output  for  the  same  input  parameters  as  Tables  3-1  and  3*2  of 

TR-741  respectively,  except  for  the  dimensionless  parameters  \,  k„,  and 

8 

k  .  A  modified  tropical  atmosphere  was  again  used.  The  two  tables  are 
for  a  high  yield  (5  MT)  and  low  yield  (20  KT)  explosion,  respectively. 
Appendix  E  gives  a  glossary  of  printout  symbols.  All  quantities  are  in 
mks  units  except  explosion  energy,  W,  which  is  in  kilotons.  For  values 
cf  physical  constants  and  initial  conditions  used  in  numerical  solution 
ol  the  cloud  equations,  see  Appendix  C  of  TR-741.  For  a  discussion  of 
appropriate  values  of  the  dimensionless  parameters  X,  k^  and  k^,  see 
Appendix  D. 
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2.4  The  Energy  Cycle  of  the  Nuclear  Cloud 

If  the  momentum  equation  (2.3.5)  is  multiplied  by  u,  the  temperature 
equation  (2.3.6)  by  c^  and  Jie  height  equation  (2.3.8)  by  g,  then  >rith 
the  turbulent- energy  equation  (2. 3. 7 )>  they  fom  a  set  of  four  energy 
equations  whose  terms  have  the  indicated  meanings. 


Kinetic  energy,  KE 


gain 
from  H 


loss  to  PE  ^ddy-viscous  inelastic- col-  dilution 
i.e.work  done  loss  to  lision  loss  to  by 
lifting  cloud  turbulence  turbulence  entrainment 


T 

=  ~  gu  -gu 
e 


..  vT  2 

2^T  U 


2  2 
u  1  dm  u  1  dm 

2  m  dt  2  m  dt 


(2.4.1) 
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Kinetic  energy  here  means  energy  of  directed  motion  of  the  cloud,  u  /2, 
turbulent  energy,  E,  is  actually  kinetic  energy  of  randomly  moving  lumps 
of  fluid  within  the  cloud.  See  Sec.  3-l»l-  Adding  the  four  equations 
gives  an  equation  for  total  energy  per  unit  mass 


Ht  (u2/2  *  V  +  E  +  8*) 


(2-4-5) 


The  dimensionless  constants  kg  and  k^  do  not  appear  in  the  total  energy 
equation;  they  only  affect  the  rate  of  transfer  between  different  forms 
of  energy. 

The  present  cloud  model,  then,  is  based  on  conservation  of  mass 
and  energy.  The  equations  for  velocity  (or  momentum  or  kinetic  energy), 
temperature  (or  enthalpy),  turbulent  energy  and  height  can  be  considered 
merely  rules  for  energy  transformation,,  These  rules  imply  an  energy 
cycle  which  may  be  diagrammed  as  follows.  The  heavier  arrows  indicate 
the  usual  direction  of  transfer  (when  T>Te  and  u=*0) 


Kinetic 
•  Energy 


Potential 

(Gravitational) 

Energy 


Turbulent 

Energy 


Enthalpy 


The  energy  flow  through  turbulence  is  always  in  the  stme  direction; 
this  part  of  the  cycle  is  intrinsically  irreversible.  The  novelty  in 
this  model  is  the  use  of  turbulence  as  a  "delay  line"  between  kinetic 
energy  of  rise  and  heat.  Although  momentum  is  lost  by  drag  or  eddy 
viscosity,  kinetic  energy  is  not  lost  but  rather  transformed. 


Because  of  entrainment,  the  cloud  cannot  be  treated  as  a  closed 
system,  with  constant  total  energy  per  unit  mass.,  except  under  special 
atmospheric  conditions.  To  examine  change  in  total  cloud  energy  (not 
per  unit  mass)  all  mass  should  be  referred  to  the  same  datum,  even  though 
the  mass  was  entrained  at  different  altitudes  and  temperatures  (at 
different  potential  energies  and  enthalpies).  As  a  datum,  we  use  sea 
level,  z  =  0,  with  the  corresponding  ambient  temperature,  Tg  =  Te  . 

Multiplying  Eq.  (2.4.5)  by  m  and  rearranging  terms,  the  rate  of 
change  of  total  energy  relative  to  this  datum  is 


d_ 

dt 


m(u^/2 


c  (T- 
P 


T  ) 
eo 


+  E  +  gz) 


=  -  T(t 

dt  '  e 


gz/cp)- 


eo 


P 

(2.4.6) 


Now  +  6Z/Cp  is  the  potential  temperature  of  the  atmosphere  at  altitude 
z,  (the  temperature  air  at  the  given  altitude  would  have  if  compressed 
adiabatically  to  sea-level  pressure),  so  that  total  cloud  energy,  referred 
to  sea  level,  increases  with  altitude  in  proportion  to  the  excess  of 
potential  temperature  over  sea-level  temperature.  We  can  now  sketch  a 
total-energy  cycle  for  the  cloud,  referred  to  the  sea-level  *4*  cum.  (See 
Fig.  2.1) 

In  an  adiabatic  atmosphere  (one  in  which  potential  temperature  does 
not  vary  with  altitude),  we  have 


so  that,  substituting  in  Eq.  (2.4.6), 

SL  ^(u£/2  +  E  ♦  Cp(T-Te0)  *  62)]  •  o 
i.e.,  total  energy,  referred  to  sea  level,  is  constant. 
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3-  TURBULETCE  IN  13iE  NUCLEAR  CLOUD 


3.1  Theoretical  Background 

Turbulent  motion  of  a  fluid  has  by  definition  a  fluctuating,  random 
nature  in  time  and  space.  This  nature  has  suggested  Fourier  analysis 
o  velocity  components  (and  other  fluctuating  properties)  into  a  spectrum 
of  motions  of  various  length  scales,  or  "eddy  sizes".  We  shall  be  con¬ 
cerned  mostly  with  small  eddies.* 

3.1.1  The  Energy  Cascade 

The  kinetic  energy  of  a  large-scale,  intensely  turbulent  flow  is 
supposed  to  be  transm:  tied  from  the  over-all  flow  to  the  largest  eddies 
(turbulent  motions)  and  through  a  cascade  of  successively  smaller  eddies* 
to  an  approximate  minimum  eddy  size,  and  finally,  to  be  dissipated  by 
viscosity  to  heat  (molecular  motions).  This  idea  dates  back  at  least  to 
L.  F.  Richardson,  c.1920.  The  cascade  picture  suggests  that  eddies 
sufficiently  far  along  the  cascade  may  be  independent  of  large-scale 
flow  conditions,  and  leads  to  the  equilibrium  hypothesis.  (See  next 
section) 

3.1.2  The  Equilibrium  Hypothesis 

During  the  transfer  cf  energy  to  smaller  and  smaller  eddies,  the 
direction  of  motion  becomes  random  and  the  concentration  of  energy  in 
particular  eddy  sizes  (wave  lengths)  is  smoothed  out  by  transfer  to 
neighboring  sizes.  Information  on  eddy  orientation  is  lost  during 
energy  transmission,  so  to  speak.  Thus,  small-scale  turbulence  is  homo¬ 
geneous  and  isotropic  even  though  the  largest  eddies  (such  as  the  vortex 
ring  of  a  nuclear  cloud)  are  not.  The  motions  at  different  small  length 


*  'Small  eddy*  is  used.  ,  .as  a  concise  term  for  a  Fourier  component 

2 

belonging  to  a  small  length  scale,  or  large  wave  number. 
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scales  are  similar  in  the  sense  that  they  differ  only  by  scale  factors 

which  a~e  independent  of  large-scale  flow  conditions.  Eddies  of  these 

sizes  are  said  to  be  in  the  equilibrium  range.  Ibis  is  "Kolmogorov's 

hypothesis"*’;  ''the  sn»ll  scale  components  of  turbulence  are  approximately 

„2 

in  statistics!  equilibrium.  Equivalent  names  for  this  hypothesis  are: 
theory  of  "locel  similarity",  of  "universal  equilibrium",  or  of 
"similarity  of  smali  eddies".  Kolmogorov's  hypothesis  applies  only  to 
very  high  Reynolds- number  floss  (Re  >10  )  such  as  a  rising  nuclear  cloud. 

While  sufficiently  small  eddies  arc  independent  of  large-scale  flow, 
all  except  very  «mall  eddies  are  independent  of  fluid  viscosity.  This 
latter  independence  leads  to  •file  principle  of  Reynolds- number  similarity 
(Ref.  10,  Sec.  5.4): 

"Geometrically  similar  flows  are  similar  at  all  sufficiently 
high  Reynolds  numbers.  The  aspects  of  the  motion  that  are 
excluded  from  the  similarity  are  quite  simply  those  that  occur 
at  low  Reynolds  numbers,  i.e.,  those  in  which  viscous  forces 
are  comparable  with  inertial  and  pressu  e  forceB," 

This  principle  is  applied  in  the  cloud  model,  for  instance,  in  the  choice 
of  characteristic  length  involved  in  momentum  and  energy  transfer.  This 
length  is  scaled  to  cloud  size,  independent  of  fluid  viscosity.  It  is 
used  again  In  developing  a  theory  of  turbulent-diffusive  dispersion  of 
particles  from  the  cloud.  (Sec.  5)  Note  that  the  principle  of  Feynolds- 
number  similarity  id  a  states nt  about  large-scale  conditions,  while 
that  of  local  similarity  is  a  statement  about  small-scale  conditions. 

The  equilibrium  range  of  the  eddy  spectrum  is  independent  of  over¬ 
all  flow  dimensions  and  thus  of  energy  input  by  large  eddies  at  one  end 
of  the  spectrum.  Furthermore,  for  extremely  large  Reynolds  numbers^1  as 

in  nuclear  cloud  turbulence,  the  larger-eddy  part  of  this  equilibrium 

2 

range,  called  the  "inertial  sub-range"  is  also  independent  of  viscous 
*  For  translations  of  Kolmogorov's  original  papers,  see  Ref.  9. 
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dissipation  by  the  smallest  eddies  (the  "dissipation  range")  at  the 
other  end  of  the  spectrum.  The  existei  ce  of  the  "inertial  range"  or 
"sub-range"  is  sometimes"^  called  "Kolmogorov's  second  hypothesis". 
Dimensional  analyf is,  applied  to  the  inertial  and  dissipation  ranges, 
predicts  the  velocity,  energy,  etc.,  at  each  wave  length.  A.,  or  wave 
number,  ky  (reciprocal  of  wave  length).  See  Sec.  3«1*3* 

For  a  quasi- s'^ady  fluid  flew,  the  rate  of  energy  dissipation  per 
unit  mass,  £,  by  the  smallest  eddies,  is  equal  to  the  rate  of  energy 
transmission  through  the  spectrum  from  the  large  eddies,  (energy- input) 
to  the  smallest,  dissipative  eddies.  That  is,  if  A.  is  an  eddy  size 
much  smaller  than  the  maximum  flow  dimensions,  and  somewhat  larger  than 
the  minimum  eddy  size,  A.  ,  then  the  rate  of  energy  transmission  past 
A.  is  equal  to  £.  (The  size  range  of  A.  might  be  about  10Xq<  A.  (  iJ/lO, 
but  see  Sec.  3«1«^)  A  quasi-steady  fluid  flow  is  one  in  which  the 
response  time  of  the  energy  spectrum  at  a  given  small-eddy  wave  number 
is  small  compared  with  the  time  scale  of  changes  in  the  overall  flow. 

The  dissipation  rate,  £,  is  determined  by  large-scale  flow  conditions, 
although  it  i.i  the  smallest  eddies  which  actually  dissipate  kinetic 
energy  to  heat.  At  this  largest  scale  of  tie  flow,  the  Reynolds  number 
is  extremely  large  and  viscosity  is  of  negligible  importance.  Then, 
using  dimensional  analysis,  £  must  have  the  form 
£  ~  u3// 

where  u  and  l  are  a  large-scale  velocity  and  length,  respectively,  since 

these  are  the  only  large-scale  quantities  out  of  which  a  dissipation 

2  3 

rate  (power  per  unit  mass,  or  length  /time-3)  can  be  famed.  Here  and 
later,  the  symbol  ~  means  "proportional  to  and  of  the  order  of  magnitude 
of".  This  type  of  dimensional  reasoning,  which  is  used  later  in  this 

ya 

report,  is  especially  popular  with  Russian  authors  (Levich,  Landau  and 
l4\ 

Lif shitz  ) . 

Writing  E  ~ 

it  is  seen  that  P.a  has  Che  dimensions  of  kinematic  viscosity,  and  u JL 


has  the  dimensions  of  shear.  Just  as  molecular  kinematic  viscosity  is 

proportional  to  the  product  of  molecular  speed  and  mean  free  path,  an 

eddy  viscosity  can  be  defined: 

v.  ,  ~  ui 

turb 

Then  we  can  write 

e  ~ Vturb  {m/I)2 

in  analogy  to  the  dissipation  rate  in  laminar  viscous  flow,  whose  form 

/^„\2 


is 


where  y  is  a  coordinate  perpendicular  to  u. 


3.1.3  Small-Scale  Turbulence  Parameters 

The  smallest  eddy  size,  XQ,  the  "Kolmogorov  microscale",  can  depend 
only  on  (molecular)  kinematic  viscosity  and  on  £,  so  that  again  by  dimen¬ 
sional  analysis, 

~  (v3/e  )* 

(The  constant  of  proportionality  jn  these  relations  is  usually  of  order 

unity. )  Thus,  the  microscale  is  only  weakly  dependent  (to  the  power) 

on  dissipation  rate.  It  is  typically  of  order  1  mm  in  severe  atmospheric 

storms,  where  e  ~  u3/i  ~  lO^cm^/sec3  or  exceptionally,  e  ~  10^  (See  for 

6  7 

instance,  Ref.  15).  Even  if  £  ~  10  or  10  ,  as  seems  possible  in  a 
megaton  nuclear  cloud,  (see  Table  2.1)  is  not  much  less  —  say  .3  mm 

=  300  p,  since  at  the  high  temperature  and/or  low  density  of  nuclear 
clouds,  the  value  of  v  will  be  greater  than  at  standard  temperature  and 
pressure.  (See  Sec.  3.2.3).  Further  predictions  of  characteristic 
quantities  in  the  inertial  and  dissipation  rangos  are  readily  derived 
(as  in  Refs.  7  and  13)  and  are  listed  here  (Twble  3.1). 

The  most  quoted  of  these  predictions  concerns  the  form  of  the 
energy  spectrum  in  the  inertial  range.  If  v  is  the  velocity  at  length 

l/O 

scale  X  in  this  range,  then  its  form  must  be  v^  ~  (ex)  '  ,  The  corres- 

2  2/3 

ponding  kinetic  energy  per  unit  mass  is  ~  v^  ~  (ex)  7  .  Then  the  energy, 


2  ^3 

B,  in  a  spectral  interval  dX  is  £(x)dx  ~  dX  =  e2/^\  1(^3dX. 

terms  of  wave  number,  k  ~  l/x,  this  gives  "Kolmogorov1 3  law": 


In 
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TABLE  3.1  SMALL  SCAI£  TURBULENCE  P. 


TABLE  3.1  SMALL  SCALE  TURBULENCE  PARAMETERS  (coot 
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£(k)dk~Ke2/3  k"5/3dk 

JlS  ~jp 

where  K  ie  a  universal  constant.  Experimental  and  theoretical 
estimates  indicate  the  coefficient  K  is  of  order  unity. 

3.1.4  Approximations  to  the  Energy  Spectrum 

12 

It  has  been  suggested  that  for  very  high  Reynolds  numbers,  the 
inertial  energy-distribution  range,  i.e.,  the  -5/3  law,  can  be  extended 
over  the  entire  spectrum,  from  the  largest  wave  lengths,  H,  to  the  micro¬ 
scale,  \  . 


3.1.5  The  Degrees  of  Freedom  of  a  Turbulent  Flew 

Associated  with  any  eddy  size,  there  is  a  "spectral  Reynolds 

number"  Re.  =  . 

\  v 

From  Table  3*1>  substituting  e  ~  u3// , 


Re. 


^3^/3  ^A/3 

~5—~  v  hS75j 


Re/? 


a) 


4/3 


where  Re  refers  to  the  overall  flow. 


Since  for  \  \  ,  Re  ■=  1,  the 

o  \ 


3/4 


ratio  of  size  of  largest  to  smallest  eddies  is  of  order  Re 

Taking  the  largest  eddies  as  of  nearly  the  same  size  as  the  overall 

3A 

flow,  there  are  then  Re  independent  microscale  eddies,  or  degrees  of 

9/4 

freedom.  In  a  three-dimensional  flcrv^  there  are  then  Re  '  independent 
eddies,  or  degrees  of  freedom.  This  suggests  the  problems  associated 
with  numerical  computation  of  turbulent  flows.  (See  Sec.  2.2.1). 


3.2  Estimates  of  Nuclear  Cloud  Turbulence 

In  this  discussion,  reference  is  made  to  the  revised  TR-741  cloud 
model  for  numerical  estimates  of  characteristics  of  explosion  clouds. 
(See  Tables  2.1  and  2.2  of  the  present  report.  They  are  revised 
versions  of  Tables  3.1  and  3.2  of  Ref.  1,  recomputed  with  the  changes 
in  the  model  mentioned  in  Sec.  2  and  specified  in  Appendix  A.) 
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3.2.1  Scales  and  Reynolds  numbers 

As  a  large-scale  characteristic  length  for  calculating  Reynolds 
number  we  use  cloud,  diameter.  The  size  of  the  largest  eddies  (such  as 
the  vortex- ring  cross  section),  used  as  a  length  scale  for  energy  dis- 
sipation  rate,  e~u */i,  is  less  than  this  diameter,  but  we  assume  that 
l  is  a  constant  fraction  of  cloud  diameter.  For  the  corresponding 
characteristic  velocity,  u,  we  use  the  "velocity  of  turbulence",  VZE  . 

In  the  early  part  of  cloud  rise  (up  to  about  60  sec.  for  5  MT), 
the  cloud  is  much  hotter  than  the  environment.  This  leads  to  different 
values  for  "external"  and  "internal"  Reynolds  numbers. 


Re 


u&> 


ext 


Re 


int 


u L 


Here,  u  and  £are  appropriate  velocity  and  length  scales.  Since,  for 
a  perfect  gas 


P„ 


T 


and 


Srftf 


Re 


then 


ext 


Re 


int 


■ft) 
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It  is  Re^^that  is  of  interest  in  cloud  turbulence.  Although  cal¬ 
culated  turbulence  velocities  ar,1  high  —  of  the  same  order  ad  rate  of 
rise  —  the  Mach  number  of  turbulence  is  (even  for  a  5  MT  shot)  never 
more  than  l/3  because  of  the  higher  speed  of  sound  in  the  hot  cloud 
than  in  the  cold  environment.  Thus,  compressible-flow  effects  can  be 
neglected . 


3.2,2  The  Turbulent  Spectrum  of  the  Cloud 

Since  nuclear- cloud  Reynolds  numbers  are  extremely  high  (for  a 
5  MT  cloud,  Re~10“°  from  30  to  600  sec.  after  burst)  one  can  cer¬ 
tainly  assume  homogeneity  and  isotropy  of  small  eddies.  Homogeneity 
is  also  consistent  with  the  "parcel"  cloud  method.  Because  of  the 
high  Reynolds  number  we  can  also,  with  some  confidence,  extend  the 
inertial  range  to  nearly  the  entire  spectrum  from  the  microscale  to 
almost  the  largest  eddies.  (Sec  3.1»*0 
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Because  the  nuclear  cloud  is  not  in  a  steady  state,  there  is  dif¬ 
ficulty  in  defining  the  dissipation  rate,  e  .  The  rate  of  "loss"  (dilu¬ 
tion)  of  energy  per  unit  mass  due  to  entrainment  may  be  larger  than 

E  dm 

viscous  loss  (i.e.,  — -rr^e),  and  'the  flux  of  energy  from  its  source 

m  cl Xt 

in  the  large  eddies  to  the  viscous  sink  is  actually  neither  independent 
of  time  nor  of  eddy  size.  Nevertheless,  we  propose  to  neglect  any  dis¬ 
torting  effect  of  entrainment  on  the  spectrum  and  estimate  turbulence 
parameters  quasi- statically  at  each  time  of  interest  on  the  basis  of 
E  and  £  .  That  is,  we  treat  the  turbulence  as  "quasi-steady"  as  defined 
in  Sec.  3.1.2.  This  treatment  is  justified  if  the  characteristic  time  of 
change  of  E,  sayf^  I  ^fjl  ^  is  much  greater  than  the  "local  transfer  time"'*'2 


'l  I  dEl-1 

'  say\E  I  dt 


of  energy  flow  through  the  spectrum,  T  (k) 


=  k3<?(k)]-1/2. 

dE  J 


There  are  several  possible  choices  for  : 

(1)  the  total  derivative  (net  rate  of  change)  of  E.  This  cannot 
be  used  at  some  early  times,  (e.g.,  near  22  sec.  in  Table  2.1)  where 

it  is  nearly  zero  because  increase  of  E  by  generation  of  turbulence  and 
decrease  by  dissipation  and  entrainment  nearly  cancel; 

(2)  the  generation  rate,  2kg  j  u2  +  - 

6  (2E)3/2 

(3)  the  dissipation  rate,  e  =  k^  — j - 

^  l 


The  characteristic  time  (|r  calculated  using  each  of  these 

choices  for  >  is  given  in  Table  3*2  and  compared  with  the  local 

transfer  time,  for  the  Table  2,1  cloud  computation.  For  the  local 
transfer  time,  r,  we  use  a  characteristic  time  based  on  similarity 
theory  (Table  3*l)»  T  =  (ek2)  ^3=  (\2/f  )^3.  This  form  for  T  can  be 
obtained  by  substituting  Kolmogorov's  law  (Sec  3»1»3)  in  x(k)  =  ^k3£(k)j 
From  Table  3*2,  it  is  evident  that,  for  typical  inertial- range  wave- 


(l  dg\  -1 
'E  dt/  ' 


lengths,  T  is  much  smaller  than 


so  that  the  equilibrium  treat¬ 


ment  of  small-scale  turbulence  is  justified  for  the  nuclear  cloud.  We 
might  say  that  if  the  local  transfer  time  is  small  compared  with  the 
characteristic  time,  turbulent  energy  is  flowing  through  the  eddy  spec¬ 
trum  like  an  incompressible  fluid. 


<H‘  ,*  *  ’’ V 


.  M  > 

?.'»»  "'.i 


‘t 


M 

ii  i  in  H  b  aMiti  ^ liilM  i  r  i  i  . 


/a 


4.  VORTEX  MOTION  IN  THE  CLOUD 

Pictures  of  nuclear  clouds  often  show  a  vortex  ring  of  hot  incan¬ 
descent  gas  within  the  rising  cloud  (see  for  instance,  Figure  2.Cf?b  of 
Ref.  3)>  and  also,  even  at  later  times  when  the  gas  is  no  longer  incan¬ 
descent,  a  general  toroidal  circulation  of  the  cloud  boundary.  These 
phenomena  appear  more  clearly  in  air  bursts  than  surface  blasts.  Sur¬ 
face  material  in  and  around  the  cloud  may  diminish  their  visibility 
or  actually  interfere  with  their  occurrence,  say  by  mass  loading  of 
the  circulating  gas.  We  will  generally  refer  to  these  phenomena  as 
"toroidal  circulation",  restricting  "vortex  ring"  to  its  classical 
hydrodynamic  meaning.  The  reasons  for  interest  in  the  circulation  are 
(l)  that  it  may  affect  cloud  rise  and  expansion  and  especially  (2)  its 
effect  on  particle  dispersion  from  the  cloud. 

Attempts  to  calculate  circulation  and  dispersion  in  terms  of  laminar 
vortex  flow  are  con^ '  lered  in  this  section.  In  the  next  section,  we 
suggest  an  alternate  t -c  aulent- diffusion  treatment  of  dispersion. 

To  calculate  the  flow  in  and  around  the  cloud  from  first  principles 
requires  a  "local"  method  (section  2.2)  but,  as  noted,  numerical  solution 
of  the  equations  of  fluid  motion  for  the  nuclear  cloud  is  not  within 

the  state  of  the  art  of  machine  computation.  Is  it  either  possible  or 
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desirable  to  use  any  of  the  classical  theory  of  vortex  rings  (see  Lamb  , 
Chapter  VII )  to  describe  tne  atomic  cloud?  The  conditions  under  whibh 
the  classical  theory  was  developed-- incompressible  steady  flow,  homo¬ 
geneous  medium,  conservation  of  vorticity,  etc. --are  not  valid  here. 

The  cloud  starts  rising  from  rest  and  input  of  turbulent  energy  (or 
generation  of  vorticity)  continues  during  the  entire  period  of  rise. 

The  streamlines  of  a  vortex  ring  flow  are  closed,  a  condition  which 

cannot  be  reconciled  with  the  mass  entrainment  characteristic  of  atomic 
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clouds.  The  nuclear  cloudis  better  described  as  a  "thermal"  (a  mass 
of  heated,  buoyant  air)  than  as  a  vortex  ring.  Thermals 
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"mix  with  the  surroundings  as  they  advance,  and 
they  grow  bigger  thereby.  After  a  time,  therefore, 
most  of  the  fluid  of  which  a  thermal  is  composed 
was  originally  motionless  in  the  surroundings  and 
most  of  the  momentum  and  vorticity  it  possesses  will 
have  been  generated  by  the  buoyancy  forces  since  its 
birth." 

"The  difference  between  a  vortex  ring  and  a  thermal 
is  that  in  the  former  case  all  the  vorticity  and 
vertical  momentum  have  been  supplied  by  an  impulse..,, 
whereas  a  thermal  is  ideally  a  configuration  of 
vorticity  and  momentum  produced  entirely  by  the 
action  of  the  buoyancy  forces.  Bie  vorticity  anu 
vertical  momentum  produced  in  the  latter  manner 
increase  in  proportion  to  the  time  for  which  the 
buoyancy  has  operated,  and  will  therefore  ultimately 
exceed  any  impressed  upon  the  fluid  at  the  beginning. 

A  vortex  ring  composed  of  buoyant  fluid  will  therefore 
ultimately  become  a  thermal  and  the  laminar  flow  will 

break  down  so  that  mixing  with  the  fluid  into  which 

1  ft 

it  is  advancing  can  take  place. "x 

From  the  preceding  quotation,  one  might  expect  that  at  least  early 
in  cloud  rise  the  vortex-ring  treatment  could  be  applicable.  This  is 

IQ 

also  suggested  by  Turner's  '  remark  on  the  rise  of  buoyant  vortex  rings, 
that  early  in  their  rise, 

"the  ring  diameter  is  increasing  faster  than  the  ... 
material  can  spread  as  the  result  of  diffusion." 

In  our  cloud  model,  this  early  rise  would  be  the  time  when  rate  of  cloud 
rise  exceeds  average  velocity  of  turbulence;  |U|>\/2E  .  (see  assumption 
4  of  the  cloud  model,  Sec.  2.3)  But  we  find  by  computation  (Tables  2.1 
and  2.2)  that  rate  of  cloud  rise  is  nearly  always  less  than  average 
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velocity  of  turbulence:  u|=v2E  except  during  the  first  few  seconds 
of  cloud  rise,  so  that  in  Turner's  terms,  cloud  expansion  is  nearly 
always  governed  by  turbulent  diffusion. 


4.1  Late  Horizontal  Cloud  Expansion 

It  has  been  suggested  that  classical  vortex  theory  be  applied  to 
the  late  horizontal  expansion  of  megaton  clouds,  using  the  method  of 
images.  The  cloud  would  be  represented  by  a  vortex  ring  approaching 
the  tropopause  from  below  while  an  Image  ring  approached  from  above. 

The  image  ring  induces  an  expansion  of  the  approaching  cloud  vortex 
ring.  This  is  mathematically  equivalent  to  treating  the  tropopause 
as  a  rigid  boundary. 

The  physical  objections  to  this  treatment  would  be  the  same  as 
those  to  any  classical  vortex  treatment  as  given  above.  Furthermore, 
the  tropopause  is  by  no  means  impenetrable,  and  the  height  of  the 
boundary  plane  would  have  to  be  specified  as  a  function  of  yield  and 
actual  tropopause  height. 

On  the  other  hand,  (Section  2.3,  assumption  4)  the  present  model 
does  predict  the  late  expansion  of  megaton  clouds,  in  effect  as  a 
turbulent  diffusion  without  invoking  a  rigid  or  yield-dependent 
tropopause. 

4.2  The  Pasted  Vortex 

Atomic  cloud  vortex  models  proposed  to  date  consist  of  an  assumed 
classical  vortex  flow  superimposed,  or "pasted",  on  a  parcel  loethod 
cloud,  without  regard  even  to  conservation  of  energy,  much  less  to 
whether  the  resulting  "pseudo- local"  velocities  satisfy  the  time- 
averaged  Navier- Stokes  equations.  Since  the  streamlines  thus  drawn 
in  the  cloud  are  intended  to  represent  average  flow  directions  in 
the  intensely  turbulent  cloud,  the  time  of  averaging  must  be  long  com¬ 
pared  with  cloud  rise.  The  time  of  averaging  then,  should  be  just 
shorter  than  the  periods  of  the  largest  eddies,  i.s.  of  the  vortex  ring. 
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Actually,  the  local  velocities  can  represent  only  one  solution  to  the 
equation  of  motion,  that  for  steady  motion  in  an  infinite  homogeneous, 
non- viscous,  incompressible  medium;  the  vortex,  however,  is  being  super¬ 
imposed  on  radically  different  conditions. 

If  the  pasted  vortex  is  to  be  used,  it  is  necessary  to  specify 

first  the  vortex  form  and  then  its  numerical  parameters.  The  two  ex- 
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tr ernes  of  vortex  structure  are  Hill's  spherical  vortex,  in  which  the 

entire  moving  spherical  fluid  body  is  the  vortex,  end  a  thin  ring 

vortex,  which  constitutes  only  a  small  part  of  the  traveling  fluid 

volume.  The  use  of  Hill's  vortex  has  been  suggested  for  cumulus  cloud 
20 

elements.  Observations  of  nuclear  clouds  indicate  that  the  thin  ring 
is  more  appropriate  than  the  spherical  vortex,  (and  that  early  in  cloud 
rise,  the  ring  cross  section  radius  is  not  more  than  l/lO  of  the  cloud 
radius). 

Given  the  general  ring-vortex  form,  its  numerical  parameters  can 
be  obtained  either  empirically  or  theoretically.  The  flow  is  completely 
determined  by  specifying  ring  radius,  vortex-core  cross- sect  ion  radius, 
and  circulation  (line  integral  of  velocity  around  a  closed  circuit 
through  the  ring).  If  cross-section  radius  is  very  small,  i.e.,  a 
thin-ring  vortex,  then  local  velocities  outside  the  core  are  completely 
determined  by  ring  radius  and  circulation.  YLtematively,  a  number  of 
local  values  of  velocity  can  be  used  to  determine  the  vortex  parameters 
and  the  flow  structure. 

4.2.1  Vortex  Parameters  From  Cloud  Filins 

In  the  empirical  method,  vortex  parameters  are  estimated  from 

ViOti-n  pictures  of  atomic  clouds.  The  motion  of  cumulus-like  lobes 

or  irojections  on  the  rising  cloud  is  supposed  to  correspond  to  the 

locil  velocities  arotind  the  vortex  ring.  This  method  has  been  well 

21 

described  by  Norment.  He  notes  that  the  theoretical  rate  of  vortex 
translation  given  by  vortex  parameters  calculated  from  the  observed 
local  velocities  does  not  generally  agree  with  observed  rate  of  cloud 
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rise.  Presumably  vortex  parameters  must  be  somehow  scaled  or  inter¬ 
polated  between  yields  and  times  for  which  films  give  information. 

The  shot  from  which  most  of  the  film  data  has  been  obtained,  shot 

3 

Grable  of  Operation  Upshot- Knothole,  was  an  air  burst.  To  apply  the 
resulting  calculated  vortex  parameters  to  fallout  particle  trajectories, 
it  must  then  further  be  assumed  that  the  presence  in  the  cloud  of  con¬ 
siderable  amounts  of  condensed  material  from  a  surface  burst  does  not 
significantly  alter  the  fluid  circulation. 

4.2.2  Specification  of  Vortex  Parameters  by  the  Cloud  Model. 

The  vortex  ring  is  certainly  the  largest,  most  permanent  eddy  in 
the  turbulent  atomic  cloud.  In  fact,  the  vortex  ring  might  be  pictured 
as  the  flywheel  maintaining  the  cascade  of  energy  to  smaller  eddies  in 
the  turbulent  spectrum  (see,  for  instance,  Hinze,  sec.  3«5)  even 
after  the  end  of  cloud  rise,  when  no  more  energy  is  supplied  to  turbu¬ 
lence.  Hiis  energy  supply  picture  is  compatible  with  the  late  horizontal 
cloud  expansion.  Now  we  propose  that  the  transfer  of  energy  fro®,  mean 
motion  to  turbulence  begins  with  the  generation  of  the  vortex,  and  that 
the  vortex  ring  be  allocated  a  fixed  fraction  of  the  total  turbulent 
energy  density  E,  say  0.2  (Ref.  11)  or  0.25  (Ref.  2).  This  proposal 
is  suggested  by  Hinze' s  statement:11 

"In  the  fully  developed  state  it  is  not  the  largest  eddies 
that  will  have  the  maximum  kinetic  energy,  but  the  eddies 
in  a  higher- wave  number  range.  Die  range  of  the  energy 
spectrum  where  the  eddies  make  the  main  contribution 
to  the  total  kinetic  energy  of  turbulence  will  be 
called  the  range  of  the  energy- containing  eddies,... 

Though  the  more  permanent  largest  eddies  contain  much 
less  energy  than  the  energy-containing  eddies,  their 
energy  is  by  no  means  negligibly  small  and  may  still 
amount  to  as  much  as  20  per  cent  of  the  total  kinetic 
energy." 
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We  can  new  outline  a  procedure  for  specifying  vortex  parameters 
using  the  cloud  model.  Our  three  input  parameters  will  be  rate 
of  cloud  rise,  u,  vertical  cloud  radius,  r,  and  rotational  (vortex) 
energy,  0.2mE,  where  m  is  cloud  mass. 

,1? 


The  total  energy  of  a  thin-ring  vortex  flow  is 

*  2por  „  3or  n 
.17 


(4.1) 


<t  is  circulation,  o>  is  the  radius 

o 


Here  we  use  Lamb's"1  notation 

of  the  vortex  ring,  a  is  the  cross-section  (vortex  core)  radius. 

This  energy  consists,  (we  assume)  of  kinetic  energy  of  rise  of  the 
1  2 

cloud,  —  mu  ,  and  rotational  energy,  0.2mE.  Then  the  energy  equation 


becomes 


1  2 

-  mu  +  0.2mE 


The  velocity  of  translation  of  the  vortex  is 


17 


u  = 


4*05“ 

o 


(In 


Rd 3r 


-£) 


(4.2) 


(4.3) 


The  radius  of  the  cloud  is  related  to  the  radius  of  the  vortex 

17 

ring  by  the  fact  that  the  stream  function  equals  zero  at  the  boundary 
of  the  "cloud",  that  is,  at  the  boundary  of  the  fluid  which  is  carried 

17 

along  with  the  vortex  ring.  This  relation  gives  a  third  equation  (Lamb, 
sec.  l6l,  Equation  ll)  so  that  k,  a,  and  cd^  can  be  found.  The  algebra 
of  the  solution  is  involved  and  only  numerical  solutions  fo.”  n,  a,  and 
uP  can  be  obtained  explicitly.  In  any  case,  for  a  thiu  ring  the  stream 
function  can  be  given  in  terms  of  *,  uP  and  u,  without  using  a,  which 
does  not  affect  velocities  outside  the  coref 

For  purposes  of  illustration,  we  taxe  a  simplified  case.  Set 

oP  =  3/5  r  as  indicated  by  Figure  3  of  Reference  21.  Treat  t  .loud 

4  3 

as  spherical,  so  that  P  =  m/(^  «  rJ)  and  drop  the  second  term  in  paren- 
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thesis  in  Equations  (4.2)  and  (4.3)  since  ujo»  a*  Then,  dividing  the 
energy  equation  (4.2)  by  the  velocity  equation  (4.3)  and  svi  Ituting 


the  given  values  of  u>  and  p, 


o  =  1.85r 


u  /2  +  0.2E 


(4.4) 


2 

l.  t  .  nearly  the  end  of  cloud  rise,  u  /2  >  0.2E  so  that  approximately 

*  =  ru  (4.5) 

19 

This  result  can  be  obtained  directly  by  dimensional  analysis:  if  « 

depends  only  on  r  and  u,  it  must  have  the  form  of  Eq.  (4.5). 

Values  of  circulation  calculated  from  Eqs.  (4.4)  and  (4.5)  for  a 

20  KT  burst  (Table  2.2)  are  rather  lower  at  early  times  (5  to  10  sec) 

21  22 

than  comparable  values  estimated  from  cloud  films.  *  '  Now,  at  these 
times,  ambient  air  is  several  times  as  dense  as  the  cloud,  and  circular, 
tion  is  changing  rapidly,  so  that  the  use  of  classical  vortex  theory  is 
on  especially  shaky  ground.  However,  it  is  reasonable  that  if  circula¬ 
tion  is  due  to  the  drag  of  the  ambient,  air  on  the  cloud,  it  should  be 
proportional  to  the  density  ratio  pg/p  =  T/Tp,  like  the  eddy- viscous  mo¬ 
mentum  loss  rate  in  the  momentum  equation  (2. 3. 5 ) .  An  appropriate  formula 
would  then  be 


*  =  ruT/Tg  (4.6) 

No  theoretical  consistency  with  Eqs.  (4.1)  -  (4.3)  can  be  claimed  for 

Eq.  (4.6):  the  latter  is  offered  as  an  intuitive  value  for  circulation 

in  a  situation  where  attempts  to  use  formulas  from  classical  vortex 

theory  seem  particularly  remote  from  reality. 

Table  4.1  gives  values  of  «  calculated  from  Eqs.  (4.4),  (4.5),  and 

(4.6),  for  several  times  of  interest  using  Tables  2.1  and  2.2  (for  5  MT 

and  20  KT  bursts,  respectively)  to  supply  values  of  r,  u,  and  T/Te>  The 

resulting  values  of  «  are  of  the  same  order  of  magnitude  as  those  estJ- 

21  22 

mated  from  cloud  films  ’  assuming  a  ring  vortex  form.  These  values 
from  films  are  widely  scattered,  some  being  smaller  and  others  larger 
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VORTEX  CIRCULATION  CALCULATED  FROM  CLOUD  MODEL 


than  those  in  Table  4.1.  They  are  for  somewhat  different  yields  and 
times  than  those  in  the  table.  Our  conclusion  is  that  if  one  wishes  to 
superimpose  an  assumed  vortex  ring  on  a  parcel-method  cloud  model,  one 
may  as  veil  derive  vortex  parameters  iron  the  model  itself  as  from  film 
observations.  Not  only  is  laborious  film  measurement  avoided,  but,  at 
least  for  Eqs.  (4.4)  and  (4.5),  the  vortex  velocity  and  energy  we  new 
compatible  with  those  of  the  cloud. 

4.3  Use  of  Kelvin's  Theorem  to  Estimate  Cloud  Circulation 

The  equations  for  the  classical  vortex  ring  are  derived  assuming 

steady  incompressible  flow.  Under  such  circumstances,  circulation  is 

constant,  hence  it  is  Illogical  to  claim  that  the  rate  of  change  of 

circulation  can  be  obtained  by  differentiating  Eq.  (4.4)  or  (4.5).  An 

expression  for  rate  of  change  of  circulation*  -rr  ,  can  be  obtained  more 

at  23 

consistently  from  a  consequence  of  Kelvin's  equation, 

^iiere  v  is  velocity,  a  is  acceleration,  and  dx  is  a  line  element.  The 
integral  is  taken  around  a  circuit  moving  with  the  fluid,  i.e.,  con¬ 
sisting  of  the  same  fluid  particles. 

Writing  the  Euler  equation  of  motion  as 
-  l 

a  u  -  vp 
P 

the  rate  of  change  of  circulation  around  the  circuit  is 

S  =  dp 


*  The  symbol  V  is  used  here  for  circulation  to  conform  with  Ref.  23, 
ifhLle  in  Sec.  4.2.2  the  symbol  «  was  used  to  conform  viw>  Ref.  17. 
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To  calculate  cloud  circulation,  we  take  a  rectangular  circuit  in  a 
vertical  plane  through  the  cloud  axis.  Ac  indicated  in  the  following 
diagram,  the  circuit  extends  vertically  from  below  the  bottom  to  above 
the  top  of  the  cloud,  and  horizontally  from  the  axis  to  well  out  in  the 
ambient  air.  ) 


Now  dp  is  zero  along  the  horizontal  parts  of  the  circuit.  Furthermore, 

the  parte  of  the  upward  vertical  path  outside  the  cloud  are  exactly 

cancelled  by  the  corresponding  parts  of  the  downward  path  at  the  same 

altitudes  so  that  the  only  possible  non- zero  contribution  to  the  line 

integral  is  from  ■‘’he  remaining  vertical  parts,  i.e.,  between  the  top 

and  bottom  of  the  cloud,  a  distanced*  *  2r.  On  each  of  these  two 

vertical  parts,  take  an  average  value  of  specific  volume,  namely,  1/p 

inside  and  l/p  outside  the  cloud,  and  using  the  hydrostatic  law,  re- 
e 

place  dp  by  -pegdz.  On  the  upward  path,  the  change  in  p  is  then 

-p  gAz  =  -2 p  gr  and  similarly  on  the  downward  path  (pressures  inside 
6  © 

and  outside  being  equal)  so  that 

AT  .  .  |T.  „eg  .  |E  peg  ]  -  2rg  -  lj  (4.9) 

Note  that  the  cloud  need  not  be  spherical;  it  is  onay  its  vertical  ex¬ 
tent,  2r,  that  matters. 


The  foregoing  derivation  makes  no  allowance  for  decrease  of  circu¬ 
lation  due  to  (l)  eddy  viscosity  or  (2)  entrainment.  Inclusion  of  these 
two  effects  violates,  respectively,  the  (l)  ideal- fluid  and  (2)  closed- 
streamline  assumptions  of  Kelvin's  theorem.  Nevertheless,  we  offer  a 
way  to  allow  for  these  effects  on  circulation.  Substitute  in  Eq.  (4.7) 
the  value  of  acceleration  given  by  the  momentum  equation  ( 2,3*5)  and 
integrate  around  the  same  circuit  as  before.  The  entrainment  rate  is 
taken  at  zero  on  the  part  of  the  circuit  outside  the  cloud.  The  resulting 
rate  of  change  of  circulation  is 


qf 

dt 


„  17  pe  \  /  pe  \  2k2v  1  dm  1 

8rL(7"1)8“vT"1i"T  U"  =  *UJ 


(4.io) 


The  first  term  in  the  brackets  on  the  right  side  of  this  equation  is  the 
ideal- fluid  term,  corresponding  to  Eq,  (4.9).  'Hie  second  and  third  terms 
represent  the  effect  of  eddy- viscosity  and  entrainment  respectively.  In 
this  heuristic  derivation,  it  seems  a  worse  abuse  of  hydrodynamic  theory 
to  introduce  entrainment  across  streamlines  than  eddy  viscosity.  Figure 

4.1  gives  numerical  values  of  circulation  from  a  hand  integration  of 
Eq.  4.10  with  and  without  the  entrainment  term,  using  data  from  Tables 

2.1  and  2.2.  The  integration  was  in  1-second  steps  to  5  sec.,  then  5- 
second  steps  to  6 0  sec.,  and  10-second  steps  thereafter. 

It  is  seen  in  Fig.  4.1  that  circulation  attains  nearly  its  maximum 
value  early  in  cloud  rise,  at  about  the  time  of  maximum  rate  of  rise. 

The  time  of  rapid  increase  in  circulation  can  be  associated  with  the 
formation  of  a  vortex  ring  from  the  initial  spherical  cloud.  This  pic¬ 
ture  is  reinforced  if  the  entrainment  term  is  omitted  -  no  entrained 
mass  crosses  the  surface  of  the  initial  cloud  even  though  this  surface 
is  deformed  from  sphere  to  toroid. 
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No  direct  experimental  confirmation  of  these  theoretical  values 

of  circulation  is  available.  Film  measurements  can  give  only  peripheral 

velocities,  not  velocities  around  the  part  of  the  circuit  inside  the 

cloud.  To  obtain  a  "measured"  circulation,  a  specific  vortex  form  (the 

21  22 

ring  vortex)  is  assumed  '  giving  a  mathematical  relation  between  per¬ 
ipheral  velocity  and  circulation.  An  indirect  confirmation  of  a  good 
choice  of  superimposed  vortex  would  be  given  by  agreement  between  the 
correnponding  "measured"  circulation  and  the  theoretical  values  given 
by  Eq.  (4.10). 

4.4  Centrifugal  Throwout  From  the  Vortex  Flow 

The  effect  of  "toroidal  circulation"  in  the  cloud  on  particle  mo¬ 
tion  may  be  determined  by  calculating  first  the  fluid  motion  and  then 
the  particle  trajectories  through  the  circulating  fluid  (provided  the 
presence  of  particles  does  not  substantially  perturb  the  flow) .  For 
sufficiently  small  particles,  the  flow  psust  the  particles  can  be  con¬ 
sidered  viscous  (Stokes'  flow);  for  larger  particles,  empirical  ex¬ 
pressions  for  drag  force  must  be  used.  In  both  cases,  of  course,  tur¬ 
bulence  is  completely  neglected  if  the  trajectories  are  calculated  from 
a  laminar  vortex  flow.  (Stokes'  law  in  itself  is  not  incompatible  with 

turbulence,  provided  we  look  only  at  the  motion  of  a  small  particle 
24  7 

within  an  eddy.  ’  )  Even  if  small  particles  are  centrifugal  distribu¬ 
tion  effect  will  be  blurred  by  atmospheric  diffusion  long  before  they 
reach  the  ground.  It  is,  therefore,  larger  particles,  say  of  diameter 
of  order  lOOOp,  for  which  toroidal  circulation  could  significantly  affect 
surface  fallout  distribution.  For  such  particles,  empirical  expressions 

for  drag  force  must  be  used.  Calculations  of  particle  trajectories  have 

£1  22 

been  carried  out  by  Norment  * 
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5.  TURBULENT  DIFFUSION  OF  PARTICLES 

If  particles  escape  from  thenucle&r  cloud  by  centrifugal  force 
as  well  as  gravity,  due  to  "centrifugal  throwout"  from  the  toroidal 
circulation  of  the  fluid,  then  the  concentrat  on  of  these  oarticles, 
immediately  after  leaving  the  cloud,  would  vary  with  distance  from 
the  axis  of  the  rising  cloud.  Their  distribution,  projected  on  the 
ground,  would  increase  with  distance  from  the  projected  cloud  center, 
while  if  particles  escaped  under  gravity  s.lone,  the  distribution  would 
be  nearly  uniform  over  the  projected  cloud  area.  The  intense  turbulence 
of  the  cloud  suggests  treating  any  non- gravitational  dispersion  of  par¬ 
ticles  as  a  turbulent  diffusion.  How  do  particles  "really" get  out  of 
the  cloud?  Perhaps  the  question  is  meaningless.  We  propose  to  use  the 
language  of  turbulence,  as  more  appropriate  than  that  of  laminar  flow, 
to  model  non-gravitational  dispersion. 

5.1  Geometrical  Effect  of  Turbulent  Diffusion 

Suppose  turbulent  diffusion  of  particles  is  uniform  over  the 
spheroidal  cloud  surface  (a  sphere  before  the  cloud  top  has  reached  the 
tropopauce;  an  oblate  spheroid  thereafter).  Only  radial  diffusion 
causes  dispersion  of  particles  from  the  cloud.  Dispersion  of  particles 
out  through  the  upper  surface  of  the  cloud  can  be  neglected  since  par¬ 
ticles  so  dispersed  are  either  recaptured  by  the  rising  cloud  or  fall 
back  into  the  cloud  after  end  of  rise.  Therefore,  only  dispersion 
through  the  lower  surface  need  be  considered. 

If  particle  flux  density,  Iq,  is  uniform  over  the  cloud  surface, 
then  its  projection,  I,  on  a  horizontal  plane  varies  with  radial  dis¬ 
tance,  x,  for  the  projected  cloud  center  as  the  ratio  of  arc  length  on 
a  cross  section  through  the  vertical  axis,  ds,  to  projected  arc  length, 
dx.  If  y  is  the  vertical  coordinate,  then  ds  =  y[{  dx)^  +  (dy)2  so 

that  for  a  spheroid  whose  vertical  cross  section  has  the  equation 
2  2 

2—  +  X—  =  1,  we  have 
d  ,  d 

a  b 
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ds  Jl-(x2/a2)(l-b2/a^) 

,-xV 


Since  the  projected  flux  density  is  proportional  to  it  becomes 

infinite  below  the  periphery  of  the  cloud,  x  =  a,  (neglecting  diffusion 

in  the  atmosphere,  wind,  etc  )  A  more  realistic  approach  is  to  divide 

the  projection  into  say  five  concentric  rings  whose  mid- radii  are  x/a  = 

0.1,  0.3,  0.5,  0.7,  0.‘b  The  corresponding  relative  projected  flux 
ds 

densities  are,  for  a  sphere,  1.01,  1.05,  1.15,  1.40,  2.2f>  and 
for  a  2/l  ellipsoid  (such  as  the  5  MT  cloud  in  Table  w.l  at  l4o  seconds) 
1.00,  1.01,  1.04,  1.11,  1.44.  Thus,  the  geometrical  model  predicts  non- 
uniform,  radially  increasing  projected  particle  flux  density.  This  is 
the  same  effect  attributed  to  toroidal  circulation. 

The  total  projected  flux  is  equal  to  the  flux  density  integrated 
over  rings  of  area  2jt  x  dx,  and  equal  to  the  total  flux  through  the 
lower  half  of  the  cloud.  For  a  spherical  cloud  of  radius,  R  =  a  =  b, * 
the  demonstration  is  particularly  simple: 


Total  projected  flux 


( 


I(x)  •  2n  x  dx  =  I 


■/ 


ds 

7R 


•  2n  x  dx  = 


ux  =  2mR  I 


=  total  flux 


5.2  Turbulent  Diffusion  Coefficients 

Turbulent  diffusion  coefficients  as  a  function  of  particle  size, 

7 

have  been  determined  in  a  study  of  relative  motion  and  coagulation  of 
particles  in  a  turbulent  gas. 

The  derivation  of  the  diffusion  coefficierts  and  other  particle 

parameters  uses  the  similarity  theory  of  turbulence  and  was  suggested 

13 

by  a  remark  of  V.  G  Levich  “  that  for  any  particle  size  there  is  a  min¬ 
imum  eddy  size  which  just  contains  the  particle  during  a  turbulent  fluc¬ 
tuation.  (Exception:  small  particles  below  a  certain  size  are  contained 

*  The  symbol  R  is  used  for  cloud  radius  in  this  section  to  distinguish 
it  from  particle  radius,  r. 
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by  all,  even  microscale  eddies).  Mean  particle  speed  relative  to  fluid 
is  a  maximum  for  this  eddy  size,  since  turbulent  velocity  increases 
with  eddy  size  (see  Table  3.1). 

By  analogy  with  the  molecular  diffusion  coefficient  in  kinetic 
theory,  the  turbulent  diffusion  coefficient  is  defined  as  the  product 
of  the  mean  free  p*<,th  and  the  speed,  relative  to  the  fluid,  of  the 

7 

particle.  The  table  diffusion  coefficients  and  other  parameters  is 

reproduced  here  (Ta'.j.e  5.1) .  As  indicated  in  the  table,  there  are  four 

flow  regimes,  which  apply  to  particles  of  increasing  size,  reading  from 

left  to  right.  The  mean  speed  of  theparticle  relative  to  the  fluid  in 

13 

the  first,  microscale  regime  is  equivalent  to  that  given  by  Levich  . 

In  each  of  the  other  three  regimes,  the  speed,  q,  is  given  by 


where  t  is  the  dissipation  rate,  p  is  density  of  the  '"rticulate 

P 

material,  and  drag  coefficient  is  a  specified  fu.  of  Reynolds 

number,  Re  =  This  equation  for  q  was  derived  by  vich  for  the 

case  CD  =  constant,  and  derived  in  Ref.  7  without  thi  rertriction. 

In  these  three  regimes,  the  forms  C  =  2*+/Re  (Steles'  law), 

CQ  =  24  Re-4  (empirical)  and  =  constant  were  used,  respectively. 

The  relaxation  time,  T,  in  each  of  the  four  regimes  is  given  by 
rf> 

t  ~  •=- *—  where  Stokes'  law  is  used  for  C_  in  both  of  the  first  two 
CDqp  D 

regimes. 

The  mean  free  path  is  qT,  which  equals  the  minimum  eddy  size  con¬ 
taining  the  particle  for  all  but  the  smallest  particles.  The  turbulent- 

2 

diffusion  coefficient  is  q  t.  For  the  smallest  particles,  contained  by 
all  eddies,  the  turbulent- diffusion  coefficient  is  less  than  the 
molecular-diffusion  coefficient,  taken  as  nearly  equal  to  kinematic 
viscosity. 
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TABLE  5.1  RELATIVE  VELOCITY  PARAMETERS  AS  A  FUNCTION  OF  PARTICLE  SIZE 
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Note  that  with  increasing  particle  size,  i.e.,  moving  from  left  to 
right  in  Table  5.1,  kinematic  viscosity  becomes  ever  less  important,  ap¬ 
pearing  raised  to  smaller  absolute  powers  when  the  mlcroscal.es  \  ,v  ,t 

o'  o'  o 

(Table  3-l)  Qre  expressed  in  terms  of  v  and  E.  Finally,  in  the  extreme 
right  column,  v  does  not  appear  at  all,  except  in  the  Reynolds  number. 
The  expressions  in  this  column,  then,  are  pure  inertial- range  forms: 
they  apply  to  interaction  of  particles  with  eddies  so  much  larger  than 
the  microscale  that  the  interaction  parameters  depend  only  on  the  dis¬ 
sipation  rate,  for  a  given  value  of  the  ratio  of  particle  density  to 
fluid  density. 

Since  there  is  a  maximum  eddy  size  of  order  l,  somewhat  less  than 
cloud  radius  R  (Sec.  3.2.1),  there  is  also  a  maximum  particle  size 
for  which  containing  eddies  exist.  The  maximum  relative  speed  for  any 
particle  size  cannot  be  greater  than  the  large-scale  turbulent  velocity 
/2E.  Since  q  is  one  component  of  relative  speed,  then  the  maximum  speed 
is 

/3q  5  /2E 


The  relative  speed  for  large  particles  given  in  Table  5.1  is 
q  =  (erpp/p)3,  and  e  is  defined  by  Eq.  (2.3.9),  where  for  a  spherical 
cloud  £  -  R.  Substituting  theBe  values  of  q  and  E  in  the  inequality. 


It  is  found  empirically  (App.  D)  that  k^  ~  .115,  so  that  the  maximum 
contained  particle  size  is  at  most 


r  *  R 


For  still  larger  particles,  the  definition  of  a  turbulent  diffusion 
coefficient  does  not  make  sense,  since  such  particles  cannot  travel 
several  mean  free  paths  within  the  dimensions  of  the  cloud.  At  any  rate 
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\ 
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their  effective  diffusion  coefficient  cannot  be  greater  than  that 
formed  from  large-scale  speed  and  length,  vSeT  R. 


5.3  Turbulent- Diffusive  Flux  of  Particles 

The  flux  of  particles  through  a  unit  area  of  cloud  surface  equals 
the  product  of  the  diffusion  coefficient  and  the  concentration  gradient 
(as  for  any  transport  process).  Hie  gradient  has  the  units  of  concen¬ 
tration  per  unit  length.  This  length  is  taken  as  proportional  to  the 
cloud  radius  (vertical  radius  when  the  cloud  is  an  ellipsoid).  Hiis  is 
the  same  characteristic  length  used  in  the  cloudi-model  equations,  Section 
2.3.  That  is,  the  principle  of  Reyn olds -number  similarity,  (Sec.  3*1)  is 
again  used  here:  since  we  are  dealing  with  a  turbulent  process,  the 
length  is  independent  of  fluid,  (molecular)  viscosity,  or  diffusivity. 
This  principle  may  not  apply  to  the  smallest  particles  mentioned  above, 
governed  by  molecular  diffusion,  but  their  diffusion  rate,  and  falling 
rate,  are  extremely  small. 

The  flux  through  a  unit  area  of  cloud  surface,  then,  is 

I0-k4D|  (5.3.1) 

where  n  and  C  are  the  concentration  and  diffusion  coefficient  of  par¬ 
ticles  of  the  given  size,  and  k^  is  a  dimensionless  constant.  The 
value  of  k^  is  not  specified,  but  by  analogy  with  the  other  turbulent 
transfer  coefficients  involving  l  in  the  cloud  equations.  Section  2.3, 
we  estimate  that 

.10  £  k^  *  .25. 


The  flux  through  the  lower  half  of  the  cloud  is  then 

\  -  |  *4  D  J 

where  S  is  cloud  surface  area. 


I 


i 


(5.3.2) 
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5.4  Gravitational  vs  Turbulent-Diffusive  Dispersion  of  Particles 

The  total  gravitational  flux  equals  the  product  of  the  horizontal 
projected  cloud  area,  the  particle  concentration  (assumed  uniform  in 
the  cloud),  and  the  particle  falling  rate,  p. 


2 

jtR  np 


(5.4.1) 


For  a  spherical  cloud,  S  4*tR  and,  substituting  this  in  Eq.  (5.3.2) 

■  \  |  <5.4.2) 

So  the  larger  the  cloud  the  less  important  becomes  particle  dispersion 
due  to  turbulent  particle  diffusion,  relative  to  gravitational  disper¬ 
sion.  The  turbulent- diffusion  coefficient  for  a  given  particle  size  in¬ 
creases  with  turbulent- energy  dissipation  rate,  and,  therefore,  with 
yield  and  finally  with  radius,  R,  but  such  dependence  does  not  signifi¬ 
cantly  alter  the  relative  importance  of  turbulent  and  gravitational 
dispersion  established  above. 


5.4.1  Large  Particles 

Consider  now  particles  of  radius,  r,  so  large  that  their  drag  coef¬ 
ficient,  Cp,  in  gravitational  fall  is  constant.  Then 

/.  n  1  i 


(!^) 


(5.4.3) 


where  p  and  p  are  partiole  and  fluid  density,  and  (p  -p  )/p  ~  p  /p . 

P  P  P^ 

Such  particle!  are  so  large  that  they  are  also  In  a  constant- drag-  coef- 

ficient  regime  with  respect  to  intensely  turbulent  motion,  (i.e.,  at 

high  dissipation  rates).  Then  (see  Table  5*1)  the  particle  Reynolds 

number  is  ever  200  and  the  diffusion  coefficient  is 


-/  ■■  S'  • 

w  ...  *  N  *\ 

,  ^  ft..  -W-  *■ 


m  i  m  v  ^  1  , 


(5.4.4) 


4/3 

D  ~  tV3  r4/3  ^ 

Tims,  taking  1.5Cp  ~  1,  (a s  is  approximately  valid  for  such  Reynolds 
numbers)  and  substituting  Eqs.  (5.4.3)  and  (5.4.4)  in  Eq.  (5.4.2) 


k4 


1/3 


e 


(5.4.5) 


The  diffusion  coefficient,  Eq.  (5.4.4),  and  therefore  the  turbulent- 
to-gravitational  flux  ratio,  Eq.  (5-4.5),  are  only  valid  for  contained 
particles  of  sizes  up  to  r  ~  Rp/p  (See  Sec.  5.2).  Substituting  the 
maximum  diffusion  coefficient,  m  r,  and  Eq.  (5.4.3)  in  Eq.  (5.4.2), 
the  flux  ratio  for  oversized,  non- contained  particles  is  at  most  about 


Comparing  Eqs.  (5-4.5)  and  (5.4.6),  it  is  seen  that  the  flux  ratio 
increases  with  increasing  r  in  the  first  equation  and  decreases  with 
increasing  r  in  the  e  ,-cond.  It  is  therefore  a  maximum  for  the  just- 
contained  particles.  In  terms  of  toroidal  circuletion  (reverting  to 
the  language  of  laminar  flow)  it  would  be  these  just- contained  particles 
for  which  centrifugal  throwout  is  most  important.  Smaller  particles  tend 
to  follow  fluid  streamlines  and  larger  ones,  to  fall  like  stones. 


For  particles  of  the  just- contained  size,  r  ~  Rp/p^,  both  flux 
ratios,  Eqs.  (5.4.5)  and  (5.4.6),  reduce  to 


'max 


R 


1/3 

17^772 


r\ 


1/2 
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where  we  have  again  used  the  definition  of  <■. ,  Eq.  (2.3.9),  and  taken 

the  turbulent  velocity  /2E  as  approximately  equal  to  rate  of  cloud 

25 

rise,  u.  But  the  Froude  number  of  the  cloud,  or  ratio  of  inertial 
to  gravitational  forces,  is 
Fr  =  u2/gR 

and  these  forces  are  nearly  in  balance  during  the  "steady"  rise  of  the 
cloud,  i.e.,  Fr  ~  1  as  is  roughly  confirmed  by  the  first  pa*,  of  Table 
2.1.  Therefore  the  maximum  flux  ratio  for  any  particle  size,  is  order 
unity 

Fr~l  (5.4.7) 

Furthermore,  fallout  particles  as  large  as  Rp/p  ,  hardly  exist.  Taking 

-  4  P 

p/p  ~  10  ,  then  even  for  a  yield  as  low  as  20  KT,  this  size  is 

r  ~  2.5  cm,  taking  R  ~  250  m  (Table  2.2). 


5.4.2  Small  Particles 

The  largest  particles  to  whose  motion  Stokes'  law  can  be  applied 

are  of  radius  r  ~  k  (p/p  )^.  These  particles  are  at  the  upper  limit  of 

^  P  '  / 

the  "Stokes- inertial"  size  range1 . 

26 

From  page  1  of  Table  2.1,  the  viscosity  of  a„r  ,  p,  and  the  defin- 

3  i 

ition  of  kinematic  viscosity,  v  =  p/p,  the  microscale  k  =  (v  /C  )4  is 

9i 

of  order  300  microns,  and  the  limiting  size  r~  k  (p/p  )4  is  about 

o  P 

30  microns.  The  diffusion  coefficient  at  this  upper  limit  is 

d  =  v^£  (5.4.8) 


while  the  falling  rate  for  such  particles  is  given  by  Stokes'  law, 
2  r2 g  Pp 

P  -  o  ~rr  ~T 


9  v  p 


(5.4.9) 
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where  again  (p  -  p)/p  ~  p  /p.  Then  the  flux  ratio,  substituting 
l  /U  ^  P 

r~\  (p/p  )  /  ,  and  Eqs.  (5 . ^ .8)  and  (5.4.9)  in  Eq.  (5.4.2)  is 

o  p 

<£  2  Ip  1/2 

,  9k  _1  _E  (5-4.10) 

g  Rg^0  1 P 

Now  we  use  the  fact  that  (Sec.  3-1-5)  ft  Re'*1  ^  and  for  the  large- 

scale  length,  l,  take  cloud  radius,  R.  Using  also  the  definition  of 


Reynolds  number  Re  ~  uR/v, 


Eq.  (5-4.10)  becomes 


(5-4.11) 


As  in  Sec.  (5-4.1),  we  substitute  the  value  of  the  Froude  number, 
2 

u  /gR  =  Fr  1  in  Eq.  (5-4.ll)  giving  for  the  ratio  of  turbulent- 


diffusive  to  gravitational  flux  of  small  particles, 


(5-4.12) 


Now  even  for  a  low  yield  cloud  (say  20  KT)  Re  ~  10^,  while  (p  /p)  ~  10^, 
,  o/2  P 

so  that  <t>/<X>  ~  10  J/  ,  i.e.,  this  ratio  is  small.  To  be  sure,  for 
^  g 

these  small  particles  gravitational  dispersion  is  itself  quite  slow. 
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6.  RESULTS  AND  CONCLUSIONS 

6 . 1  Results 

The  HU7U1  cloud  model  has  been  revised  to  take  into  account  the 
decay  of  turbulent  energy  to  heat.  The  rate  of  this  decay  per  unit  mass, 
the  turbulent- energy  dissipation  rate.  Is  the  governing  parameter  for 
small-scale  turbulent  motions.  The  model  was  shown  to  imply  a  nuclear 
cloud  energy  cycle. 

The  problem  of  toroidal  circulation  in  the  nuclear  cloud  has  been 
discussed.  It  Is  shown  that  if  this  circulation  is  to  be  represented 
by  a  classical  vortex-ring  flow  superimposed  on  the  parcel-model  cloud, 
then  the  vortex  parameters,  in  particular  the  circulation,  may  be  esti¬ 
mated  from  the  cloud  model  by  taking  vortex  rotational  energy  as  a  given 
fraction  of  the  turbulent  energy,  consistent  with  conservation  of  energy 
and  with  turbulence  theory.  The  resulting  values  of  circulation  are  of 
the  same  order  of  magnitude  as  those  calculated  from  films  of  nuclear 
clouds.  Alternately,  cloud  circulation  as  a  function  of  time  may  be 
calciilated  directly  from  the  model  by  adapting  Kelvin's  theorem.  This 
approach  is  more  consistent  with  nuclear  cloud  conditions  since  it  does 
not  require  an  assumed  steady  state  vortex  form. 

It  is  proposed  that  the  centrifugal  thrcwout  of  particles  by  the 
vortex  flow  be  represented  by  a  turbulent-diffusive  dispersion,  using  a 
previously  derived  theory  of  turbulent  motions.  For  both  large  and  small 
particles,  the  rate  of  this  turbulent- diffusive  dispersion  is  shown  to 
be  small  compared  with  the  rate  of  gravitational  fallout. 

6.2  Conclusions 

Bie  energy  dissipation  rate  and  small-scale  turbulence  parameters 
can  be  calculated  from  the  cloud  model  using  turbulent  similarity  theory. 
The  model  and  turbulence  theory  can  be  used  to  specify  vortex  parameters, 
if  a  steady  state  vortex  form  is  to  be  superimposed  on  the  cloud.  The 
model  can  also  be  used  to  calculate  cloud  circulation  as  a  function  of 


time  without  use  of  the  superimposed  vortex. 

The  uurbulent- diffusive  flux  of  particles  from  the  rising  nuclear 
cloud  is  small  relative  to  gravitational  flux,  and  can  he  ignored  for 
practical  purposes  well  within  the  accuracy  obtainable  from  current 
treatments  of  cloud  rise  aid  expansion  and  fallout- particle  dLeposition. 


APPENDIX  A 

REVISIONS  TO  THE  TR-7U1  CLOUD  MODEL 

A.l  Dissipation  of  Turbulent  Kinetic  Energy  to  neat  and  Production  of 
Turbulent  Kinetic  Energy  by  Inelastic-collision  Entrainment 

Define  the  turbulent  energy  dissipation  rate  per  unit  mass  as 

,  (2E)3/2 

€  *  ko  - 


where  k^  is  a  dimensionless  constant 

E  is  the  turbulent  energy  per  unit  mass 
I  is  a  characteristic  length 

This  definition  is  consistent  with  turbulence  theory,  (Sec.  3<l) 
The  following  changes  are  required  in  the  equations  of  TR-7^1. 
Turbulent  kinetic  energy  density,  Equation  (3*7)  of  TR-7U1: 

(l)  Subtract  €  from  right  side  (corresponding  to 
assumption  6  of  Sec.  2.3*,  dissipation  of  tur¬ 


bulent  kinetic  energy) 
u2  1  dm _ .  ^  . 


(2)  Add  2“  m  dt"  to  right  side  (corresponding  to 
assumption  7  of  Sec  2.3.,  conservation  of  total 
kinetic  energy  in  entrainment) 

(3)  Multiply  first  term  on  right  side  by  2  (this 
corrects  an  error). 

The  revised  Eq.  (3-7)  reads 

dE  T*  a  u2v  u2  1  dm  -  1  dm 

at-2k2rreyr;at-c-,:Sat 

Temperature,  Equation  (3-to)  of  TR-71*!,  "dry'  (unsaturated)  cloud: 


Add  e/cp  to  right  side,  to  obtain 


dT  g  T* 

dt  “  •  FTO  T  *  U  ' 
p  e 


r 

J  V(I)“’ 


Cp(T) 


1  dm  e 
in  dt  c  (T) 


C  .< 

,»<t'  . ..  *  J- 

t  % iiif 

■* 


* 

HaifakiCw  -W* ' 


’•Wj,  l*i 


wet"  (unsaturated)  cloud 


Add 


Temperatur 


1  ♦ 


L  £  x 

2 

o  R  T 

P  a 


e 


dT 
dt  = 


,  Fquation  (3-4w)  of  TR-7^1 
to  right  side  to  obtain: 


!-  r*  <i  *  r?>  ”  ♦ 

(x-x  )  —  +  (T-T  ) 
e  c  e 

1  dm  e 
m  dt  c 

P.  e _ 

L  P  J 

_LL 

1+ 


L26X 


It? 


P  a 


ihe  denominator  allows  for  the  increese  in  saturation  vapor  pressure 
with  temperature;  some  of  the  dissipated  energy  is  absorbed  as  heat 
of  vaporization. 


A.  2  Cloud  Form 

A  certain  altitude  z^,  the  nominal  height  of  the  tropcrpause,  iB  de¬ 
signated  such  that  the  vertical  cloud  radius  remains  constant  after  the  top 
of  the  cloud  reaches  the  trcrpopause;  r  +  z  =  z^.  Ihe  corresponding  condi¬ 
tion  in  TR-74l  applied  to  the  center  of  the  cloud:  z  t=  z^.  Thus*  in  the 
modified  treatment,  the  cloud  is  a  sphere  until  r  +  z  «=  zrJ  and  a  horizon¬ 
tally  expanding  oblate  spheroid  with  vertical  radius  r  (z^)  thereafter. 

The  theoretical  justification  is  that  the  stability  of  the  strat¬ 
osphere  restricts  the  upward  turbulent  diffusion  of  the  cloud  top.  The 
practical  reason  considers  viscous  dissipation,  i.e.,  transfer  of  energy 
from  turbulence  to  heat  (Sec.  2.3.,  assumption  6);  such  energy  transfer 
leaves  less  turbulent  energy  available  for  late  horizontal  expansion  of 
megaton  clouds  than  in  the  TR-741  model.  The  cloud-form  modification 
compensates  for  the  effect  of  assumption  6  on  late  horizontal  expansion 
and  on  vertical  oscillations,  by  lessening  dilution  of  turbulent  energy 
by  entrainment  and  starting  exclusively  horizontal  expansion  earlier. 

The  new  computed  cloud  dimensions  and  velocities  (Table  2.1)  are  now 
^  about  the  same  as  those  in  TR-741  (Table  3*l)>  which  are  in  general 

agreement  with  atomic  cloud  observations. 
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A. 3  Characteristic  Velocity 

It  was  thought  at  one  time  that  the  mean  and  turbulent  velocities 
should  be  added  to  give  the  characteristic  velocity 

v  ■  ^u2  +  2E 

This  form  has  baen  discarded  in  favor  of  v  =  max  (|u|,2E)  which  gives 
better  agreement  with  observations  of  nuclear  cloud  heights  and  dimen¬ 
sions,  and  is  consistent  with  the  rise  of  buoyant  vortex  rings. 
(Section  4.1) 

A. 4  Environmental  Conditions 
(Appendix  C.2  of  TR-74l) 

Relative  humidity  is  now  given  by  a  table  of  values,  one  for  each 
1000-meter  layer  of  the  atmosphere  up  to  28,000  meter%  and  is  taken  as 
zero  above  this  altitude,  (in  Appendix  C.2  of  TR-741,  humidity  was 
given  by  a  quadratic  polynomial  function  of  altitude.) 


A. 5  Revision  of  Volume  Equation 

Volume,  Eqs.  (3.5D)  and  (3.5W)  of  TR-741 

To  calculate  change  in  cloud  volume,  V,  the  gas  law,  p  «=  R  q(x)T, 

V  CL 

is  used  in  the  differential  form 

dV  dT  .  dq(x)  dp 
V  "  T  qW  '  p 


Substituting 

dq  =  Kl/c  -  1)/(1  +  x)2J  dx 
and  using  the  gas  law  and  the  hydrostatic  law. 


dV  ,  v  r  fi _ „  +  1  dT  l/€  -  1  dx  1  dm] 

dt  I  RaTe*  T  dt  q(x) (l+x)^  dt  m  dtj 

This  equation  replaces  both  volume  equations  of  TR-741.  Change  in  q(x) 

is  now  taken  into  account  during  the  "wet"  as  well  as  the  "dry"  period. 

<3T 

The  substitution  of  the  values  of  ^  given  by  Eqs.  (3*4)  (see  Sec.  A.l) 


is  not  made  explicitly  as  it  is  unnecessary  for  computation  and  makes  the 
equation  more  complicated  in  appearance. 
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APPENDIX  B 

ERRATA  IN  TR-74l 

3.1  Typographical  errors 

These  errata  are  typographical  only.  They  do  not  affect  derived 
equations  or  the  computer  program. 


peg'.* 


15 

lines  10-11 

Eq. 

(2.3) 

Eq. 

(2.4) 

Eq. 

(2.4) 

16 

Eq. 

(2.2) 

Eq. 

(2.8) 

17 

Para. 2,  line  4 

(1), 

line  2 

18 

Eq. 

(2.9) 

19 

Eq. 

(2.11) 

line  5 

Eq. 

(2.12) 

27 

last  line 

30 

Eq. 

for  Cp(T) 

32 

Sec. 

3.5.1 

Eq. 

for  dV)  /V 

nr 

For  "a  unit  volume  of  cloud" 

Read  "a  volume  of  cloud" 

Multiply  left  side  and  first  term  on  right  side 

by  v1 

Multiply  left  side  and  first  term  on  right  side 
by  v1 

Multiply  right  side  by 
Delete  dt  from  right  side 
For  read  2 

For  "environment"  read  "entrainment" 

For  d*  read  d^ 


_  dm 

Forat 


dS. 

peA  di¬ 


read 


dm 

dt 


PeA 


Delete  X 
For  H  read  X 
Delete  X 

For  p  /p  read  p/p 
6  6  6 

Multiply  c  (T)  by  x 


Insert  minus  sign  before  each  dT) 

m 


Eq.  for  dV)m  Multiply  right  side  by  dz 

35  Sec.  3.7  2 

Eq.  (3.7 )  Insert  l  in  denominator  under  u  v.  For  E..  read  E. 

JV  K 
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57  lines  15-16  For  Eq>  (3*6),  Sec. (3.6) 

Read  Eq.  (3-7 ),  Sec. (3.7) 

73  App.  C.3 

Eq.  for  Replace  ^  by  (l-$0 

8l  line  2  For  "the  column"  read  "ten  column" 

B.2  Errors  affecting  the  cloud  model 

page 

19  Eq.  (2.12)  For  k2  read  2kg 

35  Eq.  (3.7)  For  kg  read  2kg 

36  Eq.  (3«7A)  For  kg  read  2kg 

The  corresponding  correction  in  the  computer  program  is 

87  Statement  555*  Multiply  right  side  by  2.0. 

For  the  revised  Eq.  (3.7 ), see  App.  A.l.  For  the  revised  computer  pro¬ 
gram,  see  App.  F. 


APPENDIX  C 

SYMBOLS  USED  IN  THE  REPORT 
C.l  A  Note  on  Notation 

t  This  report  uses  hydrodynamics,  thermodynamics  and  meteorology. 

These  fields  use  the  same  symbols  for  different  quantities;  consequently, 
any  notation  used  must  violate  some  usage.  For  example,  in  meteorology 
x  and  w  are  used  for  racios  of  water  vapor  and  condensed  water  mass  to 
dry  air  mass,  respectively.  But  in  hydrodynamics,  the  velocity  compon¬ 
ents  u,  v,  w  correspond  to  the  coordinates  x,  y,  z.  Since  z  is  the 
usual  symbol  for  the  vertical  coordinate,  as  in  dp  =  -  pgg  dz,  incon¬ 
sistency  cannot  be  avoided. 


C.2  List  of  Symbols  used  in  the  present  report 

a  horizontal  semi-axis  of  cloud,  or 
acceleration,  or 
vortex  cross-section  radius 


b  vertical  semi-axis  of  cloud 

C_  drag  coefficient 


c 


P 


specific  heat  of  gas  at  constant  pressure 


D  turbulent  diffusion  coefficient 

E  turbulent  kinetic  energy  per  unit  mass 

£(k)  turbulent  kinetic  energy  per  unit  mass  per  unit  wave  number 

e  energy  dissipation  rate  per  unit  mass 

Fr  Froude  number 

f  fraction  of  explosion  energy,  W,  contained  in  fireball  at 
start  of  rise 

g  acceleration  of  gravity 

H  enthalpy 

I  turbulent- diffusive  particle  flux  density 
I  horizontal  projection  of  Iq 

kg  empirical  constant  (in  eddy  viscosity) 

k^  empirical  constant  (in  dissipation  rate) 
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h 

k 

k* 

L 

/ 

M 

m 

n 

P 

<i(x) 

R* 

R 

a 

R 

Re 

r 

S 


T* 

t 

u 

V 

▼ 

w 

y 

z 


empirical  constant  in  turbulent  diffusion 

Boltzmann's  constant;  or 
wave  number  l/\ 

latent  heat  of  evaporation  of  water 
characteristic  large-scale  length 
mean  molecular  weight  of  air 
mass  of  cloud 

concentration  of  particles  of  a  given  size 


pressure,  or 
particle  foiling  rate 

ratio  of  virtual  to  actual  temperature. 


1  +  x£ 
1  +  x 


universal  gas  constant 

R* 

gas  constant  of  air  ■  — 

a 

weighted  mean  value  of  gas  constant  -  R 

ft 

Reynolds  number 

radius  of  cloud,  or 
particle  radius 

surface  area  of  cloud 

curvilinear  coordinate 

temperature,  or 
vortex  energy 

virtual  temperature,  i.e.,  Tq(x) 
time 

vertical  velocity,  or 
large-scale  velocity 

volume  of  cloud 

characteristic  velocity,  v  «■  max  (|u|,V^  )»  or 
turbulent  velocity 

total  explosion  energy  (kilotons) 

length  coordinate 

vertical  coordinate 


1  +  X/fc 

1  +  x 
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0  ratio  of  gas  density  to  total  density  of  cloud  -  [ 

P  circulation 

C  ratio  of  molecular  weights  of  water  and  air  -  18/29 

t  energy  dissipation  rate  per  unit  mass 

k  vortex  circulation 

X  empirical  constant  (in  entrainment  rate)  or 

eddy  size 

H  viscosity 

v  kinematic  viscosity  (jx/p^) 

p  density 

t  relaxation  time  of  accelerated  particle,  or 

local  transfer  time  of  turbulent  energy 

a>o  radial  coordinate  of  vortex  ring 

cj>  total  particle  flux 

SUBSCRIPTS 

a  air  (dry  air) 

e  ambient  (environment)  conditions 

g  gravitational 

h  horizontal  (radius  of  cloud) 

o  Kolmogorov  microscale  values;  or 

initial  values 

p  particle,  or 

constant  pressure 

t  turbulent 

X  values  for  eddies  of  size  X 
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APPENDIX  D 

DIMENSIONLESS  PARAMETERS  USED  IN  THE  CLOUD  MODEL 


The  dimensionless  parameters  X,  kg  and  k^  appearing  in  the  cloud 
model  (Sec.  2.3)  are  empirical;  although  intuitive  arguments  can  be  made 
to  derive  some  particular  values  (e.g.;  k  =  .25;  see  pp.  17-18  of  Ref.  1). 


D.l  Sensitivity  of  cloud  variables  to  the  dimens  ionless  parameters 
The  maximum  rate  of  cloud  rise  depends  mainly  on  kg. 

The  maximum  height  of  cloud  center  depends  mainly  on  k. 

The  late  horizontal  radius  of  the  cloud  (if  it  reaches  the  tropo- 
p&use)  depends  almost  entirely  on  k^. 


D.2  Suggested  numerical  values  of  the  parameters 

Values  of  the  parameters  found  to  give  good  agreement  of  the  cloud 
model  with  observed  cloud  heights  and  diameters  as  functions  of  time 
are: 

4 

k  «■-  0.2.  Tuis  is  the  value  given  by  G.  I.  Taylor  in  the  earliest 
discussion  of  nuclear  cloud  rise.  In  Taylor's  report;  however,  the 
entrainment  equation  corresponding  to  Eq.  (2.3.^)  had  the  additional 
factor  T/T#  on  the  right  side.  See  discussion  in  Sec,  2.6.1  of  Ref. 
1.  Sensitivity  of  X:  0.25  is  too  large;  0.15  is  too  small. 


kg  ■  0.1  For  a  sphere,  a  given  value  of  kg  corresponds  to  a  value 
for  the  drag  coefficient  of  CD  «=  16  kg/3.  See  Eqs.  (2.13)  and 
(2.15)  of  Ref.  1,  Thus,  kg  =  0.3  is  very  nearly  equivalent  to 
Cp  =  0.5.  Sensitivity  of  kg!  0.3l25  is  too  large;  0.(775  is  too  small. 


=  O.I75.  Sensitivity  of  k^:  0.25  is  too  large,  0.125  is  too  small. 


APPENDIX  E 

GLOSSARY  OF  COMPUTER  PRINTOUT  SYMBOLS 

Symbol  Symbol 

In  in  Comment 

Printout _ Text* _ 

Output  Symbols,  in  order  of  printout 


ST 

U 

X 

T 

R 

Z 

EK 

V 

WT 

TE 

M 

ES 

P 

PW 

ED 

EPS 


t 

u 

x 

T 

r 


z 

E 

V 


m 

e_ 


P 

iL 


v  u _ v  T* 

*2  Te* 


After  SWITCH  TO  ELLIPSE  R  refers  to  hori¬ 
zontal  radius,  vertical  radius  remaining 
fixed 


Saturation  water  vapor  pressure  at  T. 
Values  of  ES  printed  for  T>373  have  no 
physical  significance 


partial  pressure  of  water  vapor  in  cloud 


rate  of  loss  of  kinetic  energy  of  rise 
due  to  eddy  viscosity,  per  unit  mass 

turbulent  energy  dissipation  rate  per  unit 

mass 


*Some  of  these  symbols  appear  only  in  Ref.  1,  not  in  the  present  report. 
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Symbol 

Symbol 

In 

in  Comment 

Printout 

Text 

Parameter  Sy 

mbola,  in  order  of  printout 

DST 

Runge-Kutta  step  Bite  between  t  -  0  and 
t  -  1 

K2 

k2 

LAMBDA 

X 

W 

w 

F 

f 

PHI 

* 

TEO 

T 

eo 

Sea-level  temperature 

CHANGE 

Value  of  ST  after  which  Runge-Kutta  step 
else  changes  to  DST2 

D8T2 

Runge-Kutta  step  sire  when  ST  >  CHANCE 

BO,  Bl. 

B2 

Coefficients  of  quadratic  polynomial  in 

T  for  c  . 

P« 

DST1 

Runge-Kutta  step  sire  for  1  <  ST  <  CHANCE 

K 

Option  in  program.  K  «  1,  3  select  the 
"alternate"  equations,  (See  TR-74l).  K  «  1, 

2  select  the  characteristic  length  /*r, 

K  «  3f  4  select  l  »  100,  This  last  value 
of  l  is  in  a  daaay  equation  in  which  some 
other  value  of  H,  such  as  a  Mach- number 
correction,  can  be  Inserted, 

Al,  A2, 

A3,  AU 

<V<Va3'aU 

Lapse  rates  in  the  4  layers  of  the  atmosphere 

Zl,  Z2, 

Z3 

*1' *2’ *3 

Dividing  altitudes  of  the  4-layer  atmosphere 

PO 

po 

Sea  level  pressure 

ZT 

Altitude  such  that  if  R  +  Z  >  ZT  (SWITCH 

TO  ELLIPSE)  cloud  vertical  radius  remains 
constant  (nominal  height  of  tropopause) 

DO,  Dl, 

Da 

Coefficients  of  quadratic  polynomial  in 

T  for  c  . 

V" 

TF 

Freesing  point  used  to  select  value  of  L 
(after  SWITCH  TO  WET) 

i 


Symbol 

in 

Printout 

Symbol 

in 

Text 

Comment 

C3 

*3 

PRINT 

End  time  of  computation 

RK3 

If  RK3  *  1/  the  momentum  equation  contains 
the  initial  virtual  mass  factor.  If  RK3  -  0 
this  factor  is  omitted. 

RLH 

Percent  relative  humidity  in  1000-meter 
layers  of  the  atmosphere,  in  ascending 
order. 

APPENDIX  F 

COMPUTER  PROGRAM  FOR  2HE  CLOUD  MODEL 

The  follawlig  FORTRAN  n  computer  program  for  the  cloud  model  ia 
the  same  is  that  given  in  TR-7^1*  except  that 

(1)  it  incorporates  the  modifications  given  in  Appendix  A, 

(2)  the  Runge-Kutta  integrator,  SUBROUTINE  RED  ILL,  is  now  written 
in  FORTRAN  II  like  the  rest  of  the  program,  instead  of  FAP, 
as  in  TR-7^1* 

The  program  computes  cloud  center  height,  rate  of  rise,  temperature, 
turbulent  energy  dissipation  rate,  etc.,  as  functions  of  time.  It  gives 
a  complete  numerical  solution  of  the  set  of  differential  equations  given 
in  TR-741  with  the  modifications  given  in  App.  A  of  the  present  report. 

The  program  was  originally  written  by  David  Hutchinson  of  CEIR,  Inc. 
and  revised  by  J.  H.  Crawford'*  and  D.  C.  Galant**  of  this  Laboratory. 

Puncning  the  Input  Data  Cards 

All  the  input  data  required  by  the  program  is  read  from  a  set  of 
six  cards;  see  the  table  below.  The  input  data  words  on  cards  1-5  of 
an  input  card  set  are  punched  as  decimal  numbers,  seven  words  in  the 
first  consecutive  ten  column  fields  per  card.  Card  3  has  an  eighth 
data  word  which  is  punched  as  an  integer  in  column  71.  Card  6  of  an 
input  card  set  has  a  different  format  which  is  specified  in  the  table. 

The  program  expects  to  process  consecutive  .nput  data  car.,  sets. 

It  halts  computation  and  calls  EXIT  when  it  reads  a  card  1  of  an  input 
data  card  set  which  has  a  negative  decimal  number  punched  in  the  third 
•rord,  i.e.,  data  word  named  T.  When  consecutive  input  data  card  sets 
are  read  by  the  program,  only  card  1  and  those  data  words  on  the  other 
cards  which  differ  froj  the  preceding  input  set  need  be  punched.  The 

*  Now  at  Computer  Usage  Development  Corp.,  Po-lo  Alto,  California 
**  Now  at  NASA  Ames  Research  Center,  Moffett  Field,  California 
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program  tests  each  input  data  word;  and  if  it  has  not  been  punched, 
replaces  it  with  the  corresponding  data  word  of  the  previous  input  data 
set.  This  means  the  initial  input  data  card  set  and  card  1  of  each  set 
must  be  punched  completely,  and  the  program  interprets  (internally)  un 
unpunched  input  data  words  as  "negative"  zeros.  The  six  data  cards 
are  punched  as  follows: 


Cols 

1-10 

11-20 

21-30 

31- 40 

41-50 

51-60 

61-70 

71-80 

Card 

1 

u 

RK3 

T 

P0 

Z 

EK 

D8T 

Card 

2 

K2 

\ 

Cl 

C2 

W 

F 

Card 

3 

TEO 

CHANGE 

BST2 

BO 

B1 

B2 

dsh 

K  (must  be 
in  col  71) 

Card 

4 

Al 

A2 

A3 

A4 

Z1 

Z2 

Z3 

Card 

5 

ZT 

DO 

D1 

32 

TF 

K3 

PRINT 

Card  6  Columns  1-8,  11- 18,  etc.,  up  to  6l-68.  Two  digit  values  of 
percent  relative  humidity  for  1000m  layers  of  the  atmosphere 
up  to  28000m.  Zero  values  mcy  be  omitted. 


Parameters  Cl  and  C2  are  not  used 


■  •  i  .*  *  >•  $•  ■ 


•  FORTRAN  CLOUDRISE 

C  PROG.  3100,  CE  !R  »  J.  CRAWFORD  (CLOUD  RISE) 

CLOUDRISE  PROGRAM  MODIFIED  BY  D  C  GALANT  25  APRIL  1965 
C  MODIFIED  7/16/65  TO  PRINT  THE  ENERGY  DISSIPATION  RATF  DC  GALANT 
C  MODIFIED  FOR  FORTRAN  INTEGRATOR  RKGILL  26  MAY  1955 

C  THIS  IS  THE  MAIN  ROUTINE 

DIMENSION  PARI  37) .  A ( 37 ) .DVBL ( 8  ) ,VRL (8) ,RKG(8) 

COMMON  A.  DWT  »0RM»0U»DX»DT  »DV»D2»DEK»  SMALL T  »DVBL  ,TE*P,P3,P1»P2*N; 

1  .El,TE2.TE3.ED.NZT.RZT.R!O.QI  .  VBL 
COMMON  EPS 

EQUIVALENCE  ( VBL ( 1  I *WT ) * ( VBL ( 2 ) ,RM I , (VBL ( 3 ) ,U»A ( 1 ) ) » ( VBL ( A ) *X  *  A ( 2 ) 

1  ) , ( VBL (  5  )  «  T  ,  A  I  3  ) I , ( VBL (  6  )  ,  V » A  (  4 ) )  * (VBL( 7) ,Z»A( 51) *(  VBL (8 ) *  E<  •  A ( 6  > 

2  ) » ( A( 7) ,DSTI , (A (8  I »RK2) • ( A ( 9 ) ,RLI »(A( 10 > ,C1) ,( A(  11 ) ,C2>  »<  A( 12) *W) 

3  *(A( 13) *  F ) , (  A  (  1 4  ) *PHI l»(A( 15) »TF0) * ( A ( 1 6 ) *  CHANGE  I »(  A( 17) .DST2) , 

4  (  A< 18) ,B0) , ( A( 19) ,B1 ) .( A(20>  »B2) , ( A( 21 ) ,DST1 ) , ( A (22 ) ,K> • ( A (23) »A1 

5  ) .1  A (24)  ,A?) , (AI25) ,A3) ,(A(26) »A4)  ,( A (27) *Z1) , < A ( 28 ) *Z2) • ( A( 29) « 

6  Z3),(A(30),P0),(A!31)*ZT)»(A(32)*D0)»(A(33)*D1),(A(34)»D2)» 

7  (A(35),TF),(A( 36 )»C3*.(A(37), PRINT)  .  (DVBLI 1 ) * DWT ) . ( DVBL ( 2 ) *DRM) , 

8  (  DVBLI 3)  ,DU) . (DVBL(4) ,DX ) , (DVBL ( 5 ) » D T ) » ( DVBL ( 6 ) *DV) , ( DVBL ( 7 )  ,DZ ) , 

9  ( DVBL (  8 ) »  DF  K  ) 

DIMENSION  RLH  (28) 

DIMENSION  IRLH  (28) 

COMMON  RLH 

PRINT  9999 
9999  FORMAT ( 1H1 ) 

C  INPUT  THE  INITIAL  CONDITIONS  0LL3WED  BY  DST  AND  PARAMETERS 

1  READ  19,  (PAR  (I) ,  I  =  1  ,  7  ) 

19  FORMAT  (7F10.0) 

C  IF  T-PARI3!  IS  NEGATIVE  IT  !$  A  SIGNAL  TO  STOP 

I F ( PAR ( 3  J  >18,4,4 
18  PRINT  9999 

CALL  EXIT 

4  READ  1900,  (PAR  (I),  I  *  8,  29),  (PAR  (I),  i  =  31,  37) 

1900  FORMAT  ( 7Fl 0 .0 /7F 10 .0 , II / 7F 1 0 .0/7F 10. 0 ) 

DO  3  I  -  1,  37 

C  IF  FIELD  IS  BLANK  USE  LAST  PARAMETER  VALUE  INPUT 

IF! PAR (  I  1)2,33,2 
B  33  D-PARI I >*20140000C700 
I F  <  D  >  3,2,2 

2  Ain  -  PARI  I  ) 

3  CONTINUE 

RK3-X 
PO  ■  V 

C  INPUT  THE  RELATIVE  HUMIDITY.  704  AND  7094  PROGRAMS. 

READ  1901,  (PAR  (I),  I  =  1,  28) 

1901  FORMAT  I  7 ( 4F  2 • 0 ,  2X)> 

DO  1905  1*1,28 

IF  (PAS  (I))  1904,  1903,  1904 
B1903  IF  (PAR  (I)  +  201400000000)  1905,  1904,  1904 

1904  RLH  ( I )  =  PAR  ( I ) 

1905  CONTINUE 

C  PRINT  THE  PARAMETERS 

PRINT  17,  (All),  I  =  7,  37) ,  RIO 

17  FORMAT (5H1DST=F8»5,5H  K2=F10.7>?H  L AMBDA-F10. 6 » 5H  C1*F  9.7.5H 

2C2-F11.7.4H  W=E 1 0, 3 , 4H  F=Fx2.8,6H  PHI=F  7.4/ 

3  7H0  TEO-F 1 0, 6 »9H  CHANGE- F7 . 3 , 7H  DST2-F 1 0, 4 , 5H  B0-F11.2  , 


4  5H  B1=F11.4.5H  B2=E12.4,7H  DST1*F7.1.4H  K=I1/6H0  A1*F8.5. 

55H  A2-F8.5.5H  A3-F8.5.5H  A4-F8.5.5H  Z1*F8.0.5H  Z2»F8.0.5H  Z3 

6*F8  «  C  •  ‘.'H  P0.F9.0/6HO  ZT-F8.0.5H  D0*F8»1.5H  D1*F8.5,5H  D2-F10. 

78. 5H  TF-F8.1.5H  C3=F8.6.8H  PRI NT  *F6. 1 »6H  RK3-F4. 1 /1H0 > 

DO  1999  I  *  1.  28 
1999  IRLH  (I)  •  RLH  (I)  +  0.3 

PRINT  1910,  (IRLH  (I).  I  -  1,  28) 

1910  FORMAT  (7H  RLH  •  7(412.  2XM 

SMALLT  *  0. 

DWT-O. 

DRM-O. 

DU*0  « 

DX-O. 

DT-O. 

DV«0. 

OZ-O. 

DEK«0. 

TO«T 

TE«TEO 

TE1«TE0-A1*Z1 
TE2«TE1-A2*(Z2-Z1 ) 

TE3»TE2-A3*(Z3-Z2 I 
XX«9. 8/287. 

P1-PO*(TE1/TEO)**(XX'A1 ) 

P2«P 1*<TE2/TE1)**(XX/A2) 

P3*P2*(TE3/TE2)**(XX/A3) 

RMA0.PHI*F*W*4.18E12/(B0*(T0-TE0)+Bl/2.*  I YO*TO-TEO*TEO) 

2  +B2/3.  * f TO*TO*TO  -  TEO*TEO*TEC ) > 

RMWO-d.-PHI  )*F*W*4.18E12/<  DO  *(TO-TEO>+  D1 / 2 .* ( TO*TO-TEO*TEO >♦ 

2  D2/3.  »(TO*T0*T0-TE0*TE0*TE0)+2500000.) 

C  SET  INITIAL  COND.  FOR  X 

X  ■  RMWO/RMAO 

C  SET  I N It (AL  FOR  V 

P«PO*l (TE0-A1*Z  ) /TEO )** ( XX/A1 ) 
V«(RMAO+RMWO)*287.*TO*(1.*29.*X/18.)/  (P*(l.  X)» 

C  SET  INITIAL  COND.  FOR  R 

R«(3.*V/ 12. 5663706) **.333333333 
RM-RMAO+RMWO 

0I».5*RM*T*( 18.+29.*X)*( l.+XE)/(TE*(18.+29.*XE)*( 1»+X ) ) 

N«  1 
WT-O. 

RZT“"1 .0 
PRINT  99 

99  FORMAT! 2X.2H ST .4X.1HU.5X.1HX.7X.1HT . 6X . 1HR • 7X » 1HZ t6X »2HEK • 4X » 1HV, 

1  8X  »2HWT .4X.2HTE.5X.1HM.1 1X.2HES  »7X 1 1HP.5X »2HPW»5X,2HED* 

2  4X.3HEPS) 

C  PRINT  INITIAL  CONDITIONS 

GO  TO  35 
8  PU*U 

PX«X 
PT-T 
PZ-Z 
PEK-EK 
PV*V 
PRM*RM 

C  TAKE  A  RUNGE-ICUTTA  STEP 


CALL  RXGILLfVBL.  DVBl.  RKG.  DST.  8) 

SMALLT  -  SMALL T  +  OST 

C  TAKE  SMALL  STEPS  INITIALLY 

IFISMALIT  -1.0)8.87.88 

87  DST-DST1 

88  IF(RZT)888,89,89 

89  R«SQRTF(3.*V/<RZT*1 7. 5663706)  > 

FPS  ■  C3*(2.*FK)**1.5/RZT 

GO  TO  35 

888  R« ( 3. *V/1 2. 5663  706  >**.333333333 

EPS  ■  C3*(2.*EK)**1.5/R 
35  PW»P*X*29«/( 18.  +29. *X) 

ES-611.*  IT/273.  )**( -5.13  )*t'XPF  I  (25. *(T-273.)  )/T) 

11  PRINT  16. 

1  smallt.u.x.t.r.z .ek.v.wt. te.rm.es. p.pw.ed.eps 

16  FORMAT (F6.1.F6.1.F7.A.F7.1 . F8 .0 » F7 .0 . F7. 0 . 1  PE9 . 2 .0PF7 . A »F6 . 1 . 

1  1PE9.2.0PF9.0.F8.0.F7.0.F8.3.1PE9.2) 

C  N«l.  DRY  MODE  N=2»  WET  MODE  N=3.  SMALL  STEP  DRY  MODE 

IF  (N-2)  150.  1 5A •  1531 

150  IF  (ES-PW)  152.152.151 

152  SAVE-DS1 
DST*0.5 

C  RESTORE  VARIABLE  VALUES  AT  START  OF  LAST  STEP 

SMALLT  »  SMALLT-SAVE 
U«PU 
X-PX 
T-PT 
Z-PZ 
EK-PEK 
V»PV 
RM-PRM 
N«3 

C  NOW  TAKE  SMALL  STEPS  UNTIL  ES  LESS  THAN  PW 

GO  TO  8 

1531  IF  ( ES-PW) 41.41.8 

41  DST-SAVE 

N*2 

PRINT  77 

77  FORMAT ( 1 4H0 SW ITCH  TO  WET) 

GO  TO  151 

154  IF(WT+. 00000001 )  153.153,151 

153  N»1 
WT*0, 

DWT-0.0 
PRINT  66 

66  FORMAT  ( 14H0SWI TCH  TO  DRY) 

C  IF  RZT  IS  POSITIVE  WE  ARE  IN  ELLIPSOIDAL  MODE 

151  IFIRZT150, 1511. 1511 

C  SWITCH  TO  ELLIPSOIDAL  IF  Z  LARGER  THAN  ZT 

CHANGED  25  APRIL  1965  TOP  OF  CLOUD  AT  TROPCPAUSE  DCG 

50  I F ( Z  +  R  -  ZT)  1511,  51,  51 

51  RZT  *R 
PRINT  52 

52  F0RMT(24H  SWITCH  TO  ELLIPSE.  R=RH ) 

1511  IF( SMALLT  -CHANGE)  14.  15.14 

15  DST-DST2 


14  CONTINUE 

IF(ABSF(T)-10.)  1 ,20,20 

20  IFIR-1.)  1.21*21 

21  IF  ISMALLT  -10.122.210.210 

210  IF(ABSF(DU)*OST-.l 1211.  22.  22 

211  IF(ABSF(UI-.l  )  1.22,22 

22  IF  ( SMALLT  -PRINT ) 1 3  •  13,  1 

13  IFIZ-10000.#W**. 25)8, 8,1 

END (0.0, 0,0,0) 
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SUBROUTINE  RKGILL 
1  (  VBL  ,  DVBL.  RKG,  H.  NNN) 

FORTRAN  RUNGE-KUTTA-GILL  INTEGRATOR  FOR  CLOUORISE  EQUATIONS 
DC  GAI.ANT  27  MAY  1965  158  DEC  STORAGE 

DIMENSION  VBL (991*  DVBLI99),  RKGI99) 

C 

COMMENCE  THE  INTEGRATION  OF  THE  CLOUORISE  EQUATIONS 
C 

CALL  DERI V 
DO  2  J  *  1#  NNN 

VBLIJ)  ■  VBL I J )  +  • 5*H*DVBL ( J ) 

2  RtCGIJ)  »  DVBL(J) 

CALL  DERIV 
DO  4  J  -  1,  NNN 

VBLIJ)  *  VBLIJ)  ♦  • 29289322 *H» I DVBL I J )  -  RKGI J )  ) 

4  RKGI J)  *  .5B578644*DVBL< J)  +  . 121 32034*RKGI J) 

CALL  DERIV 
DO  6  J  «  1.  NNN 

VBLIJ)  ■  VBLIJ)  +  1 • 70 71 068 *H# I DVBL I JI-RKGI J) ) 

6  RKGI J )  =  3.41421356*DVBL< J)-  4.121 3203*RKGI J) 

CALL  DERIV 
DO  8  J  *  1.  NNN 

8  VBLIJ)  ■  VBLIJ)  +  • 166666 67 *H» I DVBL IJ)  -  2.*RKG(J)» 

RETURN 

END  I  0.1. 0,0.0) 
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*  FORTRAN 

SUBROUTINE  DERIV 

C  VERSION  OF  DERIV  TO  BE  USED  WITH  RKGILL 
C  DERIV  WITH  CHANGES  OF  11/1/65 

DIMENSION  PARI  37)  * A ( 37  I  .OVBL  (  8 ) .VBL (8 ) *RLG( 8) 

COMMON  A.  DWT.DRM.DU.DX.OT.OV.DZ.DEK.SMALLT.DVBL.TE.P.P3.P1.P2.N. 

1  TE1.TE2.TE3.ED.NZT.RZT.RK3.QI  #  VBL 
COMMON  EPS 

EQUIVALENCE  ( VBL 1 1 1 »WT )  .  I VBL ( 2  I *RM ) . I VBL I  3 ) *U.A  ( 1 )  ) ♦ < VBL( 4 ) .X  ,A < 2 1 

1  ) .(VBL(5).T,A<3> I ♦ I VBL <  6 ) • V* A (4) ) » I VBL ( 7 J .Z »A< 5 > » •  ( VBL < 8 » #EK.A(6> 

2  I  .( A(7) .DST) . (A(0 ) »RK2) * (A (9 ) .RL)  .  (At  10 »  #C1) .( A( 11 ) »C2) •( A( 12 > .W) 

3  .( A 1 13).F> ,(A(14) .PH! )  ,(A( 15) .TEOI . ( A( 16 ) *  CHANGE ) . I  A ( 1 7 ) . DST2 » . 

4  (  A< 181 ,BO) ♦ ( A ( 19 ) .Bl ) •(  A (20) .B2> . I AI211.DST1) . (A(22) »K).( A (23) .A1 

5  ) . ( A (24) .A2) • ( At 25 ) »A3) »(A(26) »A4)  .(A(27I »Z1) »(A(28 ) *22! » (A(29) ♦ 

6  Z3) .(A(30) ,P0). ( A (31 ) *ZT) . 1  A (32) *D0>  *IA(33).D1 ) . (A(34> *D2) » 

7  (A(35)»TF)»(A(36)»C3)*(A(37)*PR!NT)»( DVBL ( 1 ) *DWT ) •  (  DVBL( 2 ) *DRM) » 

8  ( DVBL ( 3 ) »DU ) « (DVBL <  4 ) . DX  J . ( DVBL ( 5  )  t DT ) . ( DVBL ( 6 ) *DV ) * ( DVBL (7) »DZ )  . 

9  ( DVBL ( 8 ) »DE< I 
DIMENSION  RLH  (28) 

COMMON  RLH 

DZ  ■  U 

C  COMPUTE  TE  AND  P 

XX»9. 8/287. 

IF(Z-Z1>80, 80*81 

80  TE«TE0-A1*Z 
P»P0*(TE/TE0I**(XX/A1 ) 

GO  TO  89 

81  IF(2-Z2)82.82.83 

82  TE«TE1-A2*(Z-Z1I 
P-P1*(TE/TE1 )**(XX/A2) 

GO  TO  89 

83  IF(Z-Z3)84,84,85 

84  TE«TE2-A3*(Z-Z2) 

P»P2*(TE/TE2)**(XX/A3I 
GO  TO  89 

85  TE»TE3-A4*(Z-Z3) 

P-P3*(TE/TE3)**(XX/A4) 

89  CONTINUE 

C  FIND  THE  RELATIVE  HUMIDITY.  (704  AND  70941  JC 

62  L  -  XMIN1F  (Z  /  1000.0  ♦  1.00000.  28.5) 

36  XE  =  109.98  *  RLH  <L) 

1  *(TE/273.)**(-5.13)*EXPF((25.»(TE-273.) )/TE 

7)/<P*29.) 

38  CPAI»BO*(T-TF)+Bl/2.#(T*T-TF*TEI+B2/3.*  ( T*T*T-TE*TE*TE ) 

CPW*=DO  ♦  D1*T  +  D2*T*T 
CP»(B0+B1*T+B2*  T*T  +  X*CPW)/(X  1.) 

OXE-(  l.+XE)/(l.+29.*XE/18.) 

QX«( l.+29.*X/18. )/(l.+X) 

OT-T/TE 

IF(RZT)35.70.70 

35  R  *  < 3. *V/12. 56637 06 (**.333333333 

SV-3./R 
RLL  «  R 
GO  TO  49 

70  RxSQRTF(3.»V/(RZT* 12. 5663706  I ) 

ECC»SQRTF(R#R-RZT*RZT  ) /R  +  l.OE  13 
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51 

52 
60 


621 


100 


SV-3.1415926*<2.*R*R+RZT*RZT*  LOGFf  <  1  ,+FCC  » / ( 1  »-ECC  )  I/ECCI/V 
RLL  ■  RZT 
49  J«K 

EPS  ■  C3  *  12.  *  EK )  **  1.5  /  RLL 

07  •  MAX1F  (  A8SF  <U),  SORTF  (2.0  *  EIC  )) 

60  TO  (50,50.51  .51 ) ,J 
50  IF  (RZT)  53.  54.  54 

53  R2  «  R 
GO  TO  52 

54  R2  ■  RZT 

60  TO  52 
R2-100. 

GO  TO  <60,61 .60,61 ) .J 

ORM  -SV*07*( l.+29.*X/18. )*T  *  RL  *RM*QXE/( (l.*X+WT)*TE) 

00-1.0 
GO  TO  621 

61  ORM  «  SV  *  RM  *  07*  RL 

00  »  QT*OX*OXE*<  l.*X)/(  l.  +  X+VfT) 

DU  ■  ( 9. 8* (  OT  *0X  *  (l«'fX)/(l.*X*WT)*OXF  -  1.) 

2  -(00*07  *2.*RK2/R2  +  DRM  /RM)*U>* 

3  { RM+Q I  * ( 1.-RK3) I /(RM+O! ) 

M-N 

GO  TO  (100, lOi. 1001, M 

OX  — (  l.  +  Xl*(X-XE!*ORM  / (RM*( l.  +  XF* > 

DT  — OX*OT*9.8*U/CP  *  OXE  CPAI*DRM  /(CP*RM) 

♦  EPS  /  CP 

RX2  IS  K2  ♦  RL  IS  LAMBDA,  RM  IS  M 
DV  ■  V*(DRM/RM-*-DT/T  +  9.8*U/<  287,*QXE*TE)  +  (  29./18.-1.  >*DX 
1  /<QX*< l.+X|»*2ll 

GO  TO  555 

THIS  IS  THE  WET  PART 

01  «  1.  +  X*29./18. 

IF(T-TF) 102,103.103 
CI-2.83E6 
GO  TO  104 
CL-2.5E6 

02  ■  CL*X/<  287. *T  ) 

18.*Q2/(  T*29.  1 
1.  +  02 
1.+  CL*03/CP 
CL* ( X-XE ) /CP  ♦  T-TE 

•  (-0X*0T*9.8*04*U/CP  *OXE  06*DRM  /RM1/05 
EPS  /  (  CP*05  ) 

-  01* (Q3*DT  +  9,8*X*U/(2S7»*TE)*OXE) 

DV  ■  V*(DRM/RM+DT/T+9.8*(l/( 287 «*0X E*TE > + ( 29./18.-1. )*DX 
1  /(QX*< l.+X)**2l ) 

DWT  «-(  l.+X+WT)*(WT+X-XE)*DRM  / ( ( 1 .+XE }*RM )  -OX 
555  ED-RK2*  U*U*07  *00/R2  *2.0 

DEK  ■  ED  -  (  EX  -  U*U/2.»  *  DRM/RM 
1  -  EPS 

RETURN 

END  (0,0, 0,0,0  i 
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